Regresión lineal sobre una muestra que abarca muchos órdenes de magnitud

Question

La cerveza es la ley de la química dice que la absorbancia de un líquido $A$ es proporcional a la concentración de $C$, así: $$A = kC$$ El estándar cosa a hacer es preparar una serie de soluciones con concentraciones conocidas, medir la absorbancia a formar una 'curva estándar' (una curva de calibración, básicamente), y hacer una regresión lineal simple en los datos para obtener la proporcionalidad (que luego puede ser utilizado para predecir las concentraciones de desconocido soluciones).

Una manera fácil de hacer esto es comenzar con una concentración conocida y realizar una serie de dilución, lo que llevaría a 2x dilución, 4x, 8x, 16x....etc. I. e. si usted comienza con una solución de $100\mu \mathrm g/\mathrm {mL}$ le gustaría obtener soluciones con $50\mu \mathrm g/\mathrm {mL}$, $25\mu \mathrm g/\mathrm {mL}$, $12.5\mu \mathrm g/\mathrm {mL}$ etc...

Ahora, al hacer la regresión lineal, se tiene un conjunto de datos con un montón de puntos de datos a bajas concentraciones y muy pocos en concentraciones más altas. Parece mucho más natural para representar este problema en una escala logarítmica. Mi pregunta es entonces, debo estar haciendo una regresión lineal de $A$ vs $C$ o $\log A$ vs $\log C$? Al comparar los modelos, parecen dar respuestas en el mismo orden de magnitud, pero en el orden de 30% diferente.

Answer 1

Ninguno de sus modelos ajustados puede ser la correcta para la ecuación original, y su ecuación original no puede ser un modelo correcto para las variables aleatorias que se observa.

Vamos a abordar algunos de los problemas uno a la vez.

1) la ecuación original es $A=kC$, pero, por supuesto, los datos no se observa en la población los valores (de lo contrario se necesitaría solamente un $x$ $y$ en el valor, ya $k=y/x$. Claramente este no es el modelo adecuado para los datos o usted no estaría tratando de ajuste de las regresiones. Vamos a ver cómo reescribir dentro de poco.

2) Si ambas variables son observados con error, usted necesita más especializada en técnicas de regresión lineal simple de los modelos.

3) Si $C$ se observa con ningún error (o muy bajos de error relativo $A$), que presumiblemente significa un modelo algo así como decir $E(A|C) = kC$.

Entonces la pregunta es "cómo modelar la variación alrededor de la media" - necesitamos una distribución de $A$ en cada una de las $C$, o al menos alguna información sobre la distribución.

a) Como usted sugiere, usted podría ajustar un modelo lineal directamente a los datos en la escala original:

$$\text{E}(A|C)=kC$$

Esto podría adecuado si la variabilidad es casi constante en esta escala ($\text{Var}(A|C)$ casi constante).

Tenga en cuenta que este modelo, basado en su modelo de población, no tiene término de intersección; pasa por el origen.

Alternativamente, si usted tiene algún otro modelo para la variación sobre la línea, tales como $\text{Var}(A|C)\propto C^p$, entonces usted podría ajustar una regresión ponderada.

b) Como usted menciona, otro posible modelo puede ser instalado en el registro de la escala:

$$\text{E}(\log_{10}A|\log_{10}C)=\log_{10}k+\log_{10}C$$

Esto podría adecuado si la variación es casi constante en la escala logarítmica ($\text{Var}(\log_{10}A|\log_{10}C)$ casi constante) ... la que te tienden a ver si la variación de energía (un) fueron 2 o cerca de 2.

Tenga en cuenta que este modelo tiene una intercepción, pero tiene un coeficiente de la pendiente de 1. Una manera para adaptarse a este modelo de ajuste:

$$\text{E}(\log_{10}A|\log_{10}C)-\log_{10}C=\log_{10}k$$

(Dicho esto, es posible que desee comprobar un modelo más general que el original, tales como la consideración de un no-de la unidad de la pendiente en (b), por ejemplo.

La ajusta en (a) y (b) el peso de los datos de manera diferente (aunque con $p=2$ en (a) los resultados serán muy cerca), por lo que no dan respuestas idénticas. En los datos difieren en aproximadamente un 26%, lo que muestra que, para un ejemplo de la elección es bastante importante.

Si usted no tiene ningún conocimiento previo de guía para la elección, el análisis de los residuos puede ser un camino para llegar a un modelo (especialmente si usted tiene otro conjunto de datos en la que se pueden basar la elección); alternativamente, usted podría hacer la $p$ un parámetro de un modelo más grande.

(Como sucede, un poco de análisis de los residuos sugiere a mí no pueden ser posibles problemas con ambos (a) con varianza constante y con el modelo en (b); ni modelo, ni siquiera, tal vez, el más general de los modelos (a)-with0intercept y (b)-con-no-la unidad de la pendiente son necesariamente adecuados (hay una sugerencia de la curvatura sobre la línea). Una de las preocupaciones es el apalancamiento en (a), que se deriva de la 'muchos órdenes de magnitud')

4) tenga en cuenta que hay numerosos otros modelos que podrían ser instaladas.

Por ejemplo, considere el $E(A^q|C) = k^qC^q$ para algunos especificado constante $q$ (el registro del modelo corresponde a $q=0$).

Answer 2

Vamos a la física (del experimento y el aparato de medición) de la guía.

En última instancia, la absorción se determina midiendo la cantidad de radiación que pasa a través del medio y las medidas llegado a contar fotones. Cuando el medio es macroscópica, la termodinámica de las fluctuaciones en la concentración son insignificantes por lo que la principal fuente de error se encuentra en el recuento. Este error (o "ruido de disparo") tiene una distribución de Poisson. Esto implica que el error es relativamente grande en altas concentraciones poco cuando la radiación pasa a través de.

Con el suficiente cuidado en el laboratorio, las concentraciones normalmente se miden con gran precisión, así que no tendrá que preocuparse por los errores en las concentraciones.

La absorbancia en sí está directamente relacionada con el logaritmo de la medición de la radiación. Tomando el logaritmo nivela la cantidad de error en toda la posible gama de concentraciones. Sólo por esta razón, es mejor analizar la absorbancia en términos de sus valores usuales en lugar de re-expresión de los mismos. En particular, se debe evitar tomar registros de la absorción, aunque eso sería simplificar la expresión de la ley de Lambert-Beer.

También debemos estar alerta ante la posible falta de linearities. La derivación de la Ley de Lambert-Beer sugiere la absorbancia vs concentración de la curva se convierte en no lineal en altas concentraciones. De alguna forma, para detectar o prueba de que esto es necesario.

Estas consideraciones sugieren un simple procedimiento para analizar una serie de $(C_i, A_i)$ pares de concentraciones y se mide la absorbancia a:

• Estimar el coeficiente de $\kappa$ como la media aritmética de $A/C$, $\hat{\kappa} = \sum_i \frac{A_i}{C_i}$.

• Predecir la absorbancia a cada concentración en términos del coeficiente estimado: $\hat{A}(C) = \hat{\kappa}C.$

• Compruebe el aditivo de residuos $A_i - \hat{A_i}$ no lineal de las tendencias en $C_i$.

Por supuesto, todo esto es teórico y algo especulativo--no tenemos ninguna real de los datos a analizar, pero es un lugar razonable para empezar. Si se repite la experiencia de laboratorio sugiere que los datos salen de la estadística comportamientos descritos aquí, luego de algunas modificaciones de estos procedimientos sería llamado para.

---

Para ilustrar estas ideas, he creado una simulación que implementa los aspectos clave de la medición, incluyendo el ruido de Poisson y posiblemente respuestas no lineales. Ejecutando muchas veces, podemos observar el tipo de variación que es probable que se encuentre en el laboratorio. Aquí están los resultados de una simulación. (Otras simulaciones en las que puede llevarse a cabo simplemente por el cambio de la partida de semillas en el código de abajo y la modificación de diversos parámetros como desee.)

Este simulador del experimento midió la absorbancia a concentraciones de $1$$1/32$. La vertical en los diferenciales de los valores de manifiesto en el diagrama de dispersión muestran los efectos de (a) ruido de disparo en la transmisión de las mediciones y (b) ruido de disparo en la transmisión inicial de la medición de concentración cero. (Observe cómo esta realidad crea algunos negativos valores de absorbancia.) A pesar de los errores resultantes no van a tener exactamente el mismo distribuciones en cada concentración, aproximadamente igual diferenciales son evidencia empírica de que las distribuciones son lo suficientemente cerca como para ser la misma que no tenemos que preocuparnos por eso. En otras palabras, no es necesario el peso de la absorbancia de acuerdo a las concentraciones.

La línea diagonal roja ha sido estimado a partir de los 50 simulaciones. Tiene una pendiente de $\hat{\kappa}=2.13$, que difiere ligeramente de la física correcta pendiente de $2$ que fue utilizado en las simulaciones. Esta desviación es tan grande, porque supuse que no era muy poca radiación a medir; el máximo de fotones en el recuento sólo $1000$. En la práctica, el máximo de la cuenta puede ser de muchos órdenes de magnitud mayor que este, que conduce a la alta precisión de la pendiente de estimaciones -, pero entonces no podríamos aprender mucho de esta figura!

El histograma de los residuos no se ven bien: es sesgada a la derecha. Esto indica algún tipo de problema. Que el problema no viene de la asimetría en los residuos en cada concentración; más bien, se trata de una falta de ajuste. Que es evidente en la boxplots a la derecha: a pesar de que los primeros cinco de ellos de la línea casi horizontal, la última de ellas--en la concentración más alta--difiere claramente en la ubicación (es demasiado alta) y la escala (es demasiado largo). Este es el resultado de una relación no lineal de la respuesta que he construido en la simulación. A pesar de la falta de linealidad está presente a lo largo de toda la gama de concentraciones, tiene un efecto apreciable en la concentración más alta. Esto es más o menos lo que iba a suceder en el laboratorio, también. Sin embargo, con sólo una calibración disponible no pudimos sacar esa boxplots. Considerar el análisis de múltiples independiente que se ejecuta si no linealidad puede ser un problema.

---

La simulación se realizó en R. Los cálculos con los datos reales, sin embargo, son sencillos de llevar a cabo a mano o con una hoja de cálculo: sólo asegúrese de comprobar los residuos de la no linealidad.

#
# Simulate instrument responses:
#   `concentration` is an array of concentrations to use.
#   `kappa` is the Beer-Lambert law coefficient.
#   `n.0`   is the largest  expected photon count (at 0 concentration).
#   `start` is a tiny positive value used to avoid logs of zero.
#   `beta`  is the amount of nonlinearity (it is a quadratic perturbation
#           of the Beer-Lambert law).
# The return value is a parallel array of measured absorbances; it is subject
# to random fluctuations.
#
observe <- function(concentration, kappa=1, n.0=10^3, start=1/6, beta=0.2) {
  transmission <- exp(-kappa * concentration - beta * concentration^2)
  transmission.observed <- start + rpois(length(transmission), transmission * n.0)
  absorbance <- -log(transmission.observed / rpois(1, n.0))
  return(absorbance)
}
#
# Perform a set of simulations.
#
concentration <- 2^(-(0:5)) # Concentrations to use
n.iter <- 50                # Number of iterations
set.seed(17)                # Make the results reproducible
absorbance <- replicate(n.iter, observe(concentration, kappa=2))
#
# Put the results into a data frame for further analysis.
#
a.df <- data.frame(absorbance = as.vector(absorbance))
a.df$concentration <- concentration # ($ interferes with TeX processing on this site)
#
# Create the figures.
#
par(mfrow=c(1,3))
#
# Set up a region for the scatterplot.
#
plot(c(min(concentration), max(concentration)), 
     c(min(absorbance), max(absorbance)), type="n",
     xlab="Concentration", ylab="Absorbance",
     main=paste("Scatterplot of", n.iter, "iterations"))
#
# Make the scatterplot.
#
invisible(apply(absorbance, 2, 
                function(a) points(concentration, a, col="#40404080")))
slope <- mean(a.df$absorbance / a.df$concentration)
abline(c(0, slope), col="Red")
#
# Show the residuals.
#
a.df$residuals <- a.df$absorbance - slope * a.df$concentration # $
hist(a.df$residuals, main="Histogram of residuals", xlab="Absorbance Difference") # $
#
# Study the residual distribution vs. concentration.
#
boxplot(a.df$residuals ~ a.df$concentration, main="Residual distributions",
        xlab="Concentration")