Descargo de responsabilidad
Soy nuevo en este sitio, relativamente nuevo en R (dos semanas de aprendizaje), tengo sólo un conocimiento muy básico en la estadística así que lo siento si estoy haciendo un error tonto allí o hacer mala pregunta o algo.
Tampoco sé cómo incrustar bien mi conjunto de datos en el post (he buscado en meta sin encontrar nada al respecto), así que he compartido un enlace desde Google Drive (si conoces una forma mejor, por favor, házmelo saber y puedo cambiarlo o hacerlo tú mismo si tienes derecho a editar).
Explicación
En la asignatura de Laboratorio físico que estoy cursando he realizado un experimento para determinar el coeficiente de atenuación de un medio óptico utilizando filtros ópticos de diferentes grosores (con este espectrómetro http://www.vernier.com/products/sensors/spectrometers/visible-range/v-spec/ ).
Generó datos que puedes descargar aquí para reproducir mis resultados: https://docs.google.com/uc?authuser=0&id=0B5x0REqWCRBuaHlXc2JWZ19sWE0&export=download
Mi objetivo es determinar $\kappa(\lambda)$ , donde $\kappa$ es el coeficiente de atenuación y $\lambda$ de onda, utilizando la ley de Beer-Lambert que dice ( $l$ es el grosor): $$\theta=x_0 \cdot \exp(-\kappa \cdot l)$$
Tengo 5 "curvas" para diferentes espesores de filtros ópticos (de 1 a 5 mm) y ajusto el modelo exponencial para cada longitud de onda $\lambda$ que me da la $\kappa$ para esa longitud de onda.
He probado esto:
Modelo linealizado
Tomando el logaritmo de la ecuación anterior y alguna manipulación algebraica obtenemos $$ \ln \left( \frac{x_0}{\theta} \right) = \kappa \cdot l$$
representado por el comando R lm( I(log(x0/temp.theta)) ~ l + 0 )
Modelo no lineal
Tomo la ley de Beer-Lambert (primera ecuación) y la modifico por nls(temp.theta ~ I(x0 * exp(-k*l)) + 0 , start= list(k = k.l)) donde se obtiene del modelo linealizado k.l = coef(lm(...)) .
Código de reproducción
He puesto x0 = 100 desde $\theta$ se mide en porcentajes.
x0 = 100 #x=0 intercept of the exp function
data = read.delim2("export2.csv")
l = data$l
l = l[-which(is.na(l))]
data$l = NULL
data$kappa.l = NA
data$sd.l = NA
data$kappa.nl = NA
data$sd.nl = NA
for(i in 1:nrow(data)){
temp.theta = as.numeric(data[i,-which(names(data) %in% c("lambda","kappa.l","kappa.nl","sd.l","sd.nl"))])
temp.lm = lm( I(log(x0/temp.theta)) ~ l + 0 )
k.l = coef(temp.lm)
data[i , "kappa.l"] = k.l
data[i , "sd.l"] = coef(summary(temp.lm))[ ,"Std. Error"]
temp.nls = nls(temp.theta ~ I(x0 * exp(-k*l)) + 0, start = list(k = k.l))
data[i , "kappa.nl"] = coef(temp.nls)
data[i , "sd.nl"] = coef(summary(temp.nls))[ ,"Std. Error"]
}Preguntas
Cuando visualizo el resultado obtengo el siguiente gráfico. Surgen estas preguntas:
- ¿Por qué los ajustes son tan diferentes en torno a la longitud de onda de 400 nm?
- ¿Cuál se ajusta mejor a los datos?
- ¿O es que los datos están tan lejos del comportamiento exponencial ahí que tengo que cortarlos donde se encuentran?
- ¿Cuál de estos modelos es estadísticamente (más) correcto?
- ¿Puedo mejorar los ajustes añadiendo pesos en función de la incertidumbre de la medición (que creo que es del 5 % según se indica en la web de Vernier)?
el ajuste con valores más altos es el no lineal 
