Resumen
En este documento se presenta el reporte reproducible de la evaluación estadística del proceso de adaptación de una metodología de cuantificación de silicio en extractos de fertilizantes.
La extracción es asistida por microondas en medio alcalino y las mediciones se realizan mediante espectrofotometría.
Resultados
Primera aproximación
###Construcción de la curva
Para los niveles de la curva de calibración se llevaron a cabo mediciones de absorbancia en \(\lambda_{max} = 660\; nm\). Cada medición se realizó por triplicado con tres volúmenes de celda diferentes provenientes de la misma solución de calibración.
C.curva <- c(0, 1,2,5,10,20,50,75,100) # Las concentraciones preparadas
Abs.curva <- c(0.000,0.013,0.026,0.075, mean(0.179,0.123,0.140), 0.269, mean(0.623,0.628,0.625),mean(1.107,1.112,1.105),mean(1.215,1.213,1.219)) # Las absorbancias leídas, los valores que no son el resultado de un promedio tuvieron 3 réplicas iguales
curva <- data.frame(Si.mg.L=C.curva, Absorbancia.u.a.=Abs.curva)
print(curva)## Si.mg.L Absorbancia.u.a.
## 1 0 0.000
## 2 1 0.013
## 3 2 0.026
## 4 5 0.075
## 5 10 0.179
## 6 20 0.269
## 7 50 0.623
## 8 75 1.107
## 9 100 1.215A partir de los resultados experimentales, es posible crear un modelo
de regresión lineal simple mediante la el uso de la función
lm() :
modelo <- lm( Absorbancia.u.a. ~ Si.mg.L , data = curva) # Se crea el modelo lineal mediante OLS (Ordinary Least Squares Regression)Una vez creado el modelo, podemos revisar la distribución de los
residuales y los parámetros de regresión mediante el uso de la función
summary():
summary(modelo)##
## Call:
## lm(formula = Absorbancia.u.a. ~ Si.mg.L, data = curva)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.083154 -0.014577 -0.014249 -0.002292 0.129740
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.0145770 0.0268986 0.542 0.605
## Si.mg.L 0.0128358 0.0005908 21.725 1.1e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.06188 on 7 degrees of freedom
## Multiple R-squared: 0.9854, Adjusted R-squared: 0.9833
## F-statistic: 472 on 1 and 7 DF, p-value: 1.104e-07En la pantalla anterior se muestran los parámetros del modelo de regresión entre la absorbancia y la concentración de las muestras.
También es posible graficar la línea de regresión junto con los datos experimentales:
library(ggplot2)
g <- ggplot(curva, aes(Si.mg.L, Absorbancia.u.a.))
g + geom_point() + geom_smooth(method ="lm")## `geom_smooth()` using formula 'y ~ x'
Análisis de residuales para confirmar supuesto de homocedasticidad
El primer elemento es un resumen de la distribución de los
residuales, los cuales para cumplir el supuesto de homocedasticidad
deben tener varianza constante y promedio cero. Al imprimir el objeto
modelo, la distribución de los residuales se presenta como
cuartiles, siendo el segundo cuartil la mediana, un estimador robusto
del valor central de una ditribución. Es posible ver los cuartiles de
una manera más gráfica mediante un gráfico boxplot:
residuales <- modelo$residuals
boxplot(residuales)
abline(h = median(residuales), col = "red", lty = 2)
median(residuales)## [1] -0.01424855En el gráfico boxplot, se crea una caja cuyo límite inferior es el primer cuartil, es decir que dentro de esta distribución de resultados experimentales el 25% de los datos son menores al valor que corresponde al primer cuartil (-0.014577). Bajo la misma lógica se pueden evaluar los valores para la mediana (-0.014249) y el tercer cuartil (-0.002292).
La distancia entre el tercer y el primer cuartil, es llamada el rango inter-cuartil (QIR)
\[\begin{align} IQR \; = \; Q_3 - Q_1 \end{align}\]
Los bigotes son dibujados de la forma:
el bigote inferior va desde el valor de la distribución mínimo que esté entre el primer cuartil menos 1.5 veces el rago inter-cuartil \(Q_1 - 1.5IQR\) y el primer cuartil \(Q_1\).
el bigote superior va desde el valor de la distribución máximo que esté entre el tercer cuartil \(Q_3\) y el tercer cuartil mas 1.5 veces el rago inter-cuartil \(Q_3 + 1.5IQR\).
En el boxplot que corresponde a la distribución de el primer modelo de regresión, el rango intercuartil es tan pequeño que no hay valores dentro de los rangos establecidos para generar los bigotes. Esto indica que la distribución es estrecha y puede presentar curtosis, o peso en los extremos. Esto puede deberse a valores atípicos, los cuales son detectados al estar por fuera de los bigotes del boxplot.
Podemos generar un boxplot para una distribución normal para comparar cómo debería verse la distribución si fuese normal:
set.seed(2)# para que los calculos de valores aleatorios sean reproducibles cada vez que se corra el código de este informe
dist <- rnorm(9, mean=0 , sd= sd(residuales)) #
boxplot(dist)
En el gráfico boxplot para una muestra aleatoria de una distribución normal con la misma desviación presentada por los residuales experimentales es posible observar como se dibujan los bigotes.
Para confirmar la presencia de curtosis, pueden utilizarse herramientas gráficas, como los gráficos de densidad, y numéricas como el cálculo de momentos centrales para curtosis y oblicuidad.
library(ggplot2)
g <- qplot(residuales)
g + geom_density(alpha = 0.4,)## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
library(moments)
kurtosis(residuales)## [1] 4.18668Los resultados del cálculo para el momento central de curtosis toman valores entre 0 e \(\infty\). un valor mayor a 3 indica una distribución “punteaguda” (Henderson, 2005) o leptocurtica en la cual no se presentan hombros. El calculo del momento central de curtosis nos permite rechazar la hipótesis de que la dispersión de la distribución se debe a la presencia de hombros en la misma, y confirmar que se debe a datos atípicos.
Para evaluar que tanto se desvía la distribución de la normalidad, es posible realizar un gráfico cuantil-cuantil, en el cual se comparan los cuantiles para una distribución normal, con los correspondientes a la distribución empírica experimental. Si la distribución se aproxima a la normalidad, los valores de los cuantiles deberían ser los mismos para la muestra experimental y para la distribución normal, entonces, al hacer un gráfico de dispersión de los cuantiles experimentales en función de los cuantiles teóricos, los puntos deben yacer sobre sobre una línea recta.
library(car)## Loading required package: carDataqqPlot(residuales, envelope = 0.95)
## [1] 8 9De manera análoga a los resultados en el boxplot, hay una concentración de datos en un rango estrecho de la distribución, que yacen dentro del intervalo del 95% de confianza para la regresión entre los residuales experimentales y los teóricos. Los valores que yacen fuera del intervalo de confianza son valores atípicos. En este caso son los elementos 8 y 9 de la curva de calibración:
curva[c(8,9),]## Si.mg.L Absorbancia.u.a.
## 8 75 1.107
## 9 100 1.215Estos niveles son los únicos niveles de la curva que tienen una absorbancia mayor a 1, y están por fuera del rango lineal de la absorbancia (Mayerhöfer 2020).
Pruebas de hipótesis para los estimadores de los parámetros de regresión
En la impresión del contenido del objeto modelo, R nos
brinda los resultados de las pruebas de hipótesis para una regresión
Prueba t individual para los coeficientes intercepto y pendiente
En la sección coefficients de la impresión del contenido
del objeto lm se presentan los coeficientes de regresión
para el intercepto y la pendiente. La fila titualda con el valor de x,
en este caso Si.mg.L corresponde a la pendiente.
summary(modelo)##
## Call:
## lm(formula = Absorbancia.u.a. ~ Si.mg.L, data = curva)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.083154 -0.014577 -0.014249 -0.002292 0.129740
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.0145770 0.0268986 0.542 0.605
## Si.mg.L 0.0128358 0.0005908 21.725 1.1e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.06188 on 7 degrees of freedom
## Multiple R-squared: 0.9854, Adjusted R-squared: 0.9833
## F-statistic: 472 on 1 and 7 DF, p-value: 1.104e-07Teniendo como ecuación de regresión:
\[\begin{align} \hat{y_i} = b_0 + bx_i \end{align}\]
Y los residuales para cada nivel \(i\) como:
\[\begin{align} e_i = y_i - \hat{y_i} \end{align}\]
La aproximación de regresión simple, la cual minimiza la suma de los residuales al cuadrado o Ordinary Least Squares Regression (OLS) tiene la ventaja de que podemos obtener una ecuación para el resultado de cada parámetro de regresión.
Para la pendiente:
\[\begin{align} b = \frac{\sum_{i=1}^n (x_i - \bar{x})(y_i - \bar{y})}{\sum_{i=1}^n (x_i - \bar{x})^2} \end{align}\]
y para el intercepto:
\[\begin{align} b_0 = \bar{y} - b \cdot \bar{x} \end{align}\]
De este modo se obtienen los parámetros \(b_0 = 0.0133425\) y \(b = 0.0124046\). Con cada parámetro está asociado un error estándar. Este error estándar se utiliza para transformar el valor real a la escala de una distribución normal estándar, en la cual se hace la prueba de hipótesis. Para hallar el error estándar de cada parámetro de regresión es necesario obtener el error estandar de los residuales, una medida de la dispersión de los valores predichos respecto a los experimentales:
\[\begin{align}
S^2_e = \frac{1}{n-2} \sum_{i=1}^n ( y_i - \hat{y_i})^2 = \frac{1}{n-2}
\sum_{i=1}^n e^2
\end{align}\]
En donde los grados de libertad son n-2 ya que estimamos dos parámetros
del modelo. Para definir un intervalo de confianza en el nivel p para
los parámetros, se necesitan los valores de t y el error estandar de
cada parámetro:
\[\begin{align} b_0 \pm t_{n-2;p} \cdot s(b_0) \end{align}\]
\[\begin{align} b \pm t_{n-2;p} \cdot s(b) \end{align}\]
Con los errores estándar de \(b_0\) y \(b\):
\[\begin{align} s(b_0) = s_e \cdot \sqrt{ \frac{\sum_{i=1}^n x_i^2}{ n \sum_{i=1}^n (x_i - \bar{x})^2 } } \end{align}\]
\[\begin{align} s(b) = \frac{s_e}{\sqrt{\sum_{i=1}^n (x_i - \bar{x})^2}} \end{align}\]
Ya que R nos brinda al imprimir el contenido del modelo con
summary(modelo) el resultado del cálculo de los errores
estándar en la columna Std. Error y el valor de t calculado
para la prueba de hipótesis en la columna t value podemos
calcular los intervalos de confianza con un nivel de significancia de
0.05.
Q <- qnorm(p=1-0.05/2)
intervalo_intercepto <- Q *0.0268986
intervalo_pendiente <- Q*0.0005908
print(round(c(modelo$coefficients[1],modelo$coefficients[2]),3))## (Intercept) Si.mg.L
## 0.015 0.013print(round(c(intervalo_intercepto,intervalo_pendiente),3))## [1] 0.053 0.001De esta manera podemos reportar los parametros de regresión con intervalos de confianza : \(b_0 = 0.015 \pm 0.053\) , \(b = 0.013 \pm 0.001\). Para reportar de una manera correcta los parámetros del modelo de regresión, hace falta reportar el resultado de la prueba de hipótesis para la diferencia significativa de cero:
Intercepto: La probabilidad de que una diferencia de cero igual o superior a la obtenida en este exprimento se vuelva a presentar por casualidad, es \(p= 0.605\), o una probabilidad del 60%. Lo cual indica que el intercepto no es significativamente diferente de cero y no es necesario para el análisis de regresión.
Pendiente: La probabilidad de que una diferencia de cero igual o superior a la obtenida en este exprimento se vuelva a presentar por casualidad, es \(p= 1.1\cdot10^{-7}\), o una probabilidad del \(1.1\cdot10^{-5}\)%. Lo cual indica que la diferencia de cero se da por una relación causal entre los cambios de señal y los cambios de concentración. La pendiente Sí es significativamente diferente de cero y es necesaria para el análisis de regresión.
En este momento es posible presentar los resultados para los parámetros de regresión de una manera adecuada: \(b_0 = 0.015 \pm 0.053 \;(p=0.605)\), \(b = 0.013 \pm 0.001 \;(p=1.1\cdot 10^{-7})\).
Representación de los datos por el modelo
La representación de los datos experimentales por el modelo de
regresión lineal se puede evaluar con el coeficiente de regresión \(R^2\) de pearson, el cual R presenta por
defecto como Multiple R-squared, ya que también se puede
calcular para modelos con varias variables independientes. En este caso,
en regresión bivariada es el mismo que el coeficiente de pearson, e
indica el porcentaje de representación de los datos usados para
construir el modelo por la ecuación obtenida. un \(R^2 = 0.9854\) indica una representación
del 98.54% de los datos experimentales:
round(9*0.98)## [1] 9Para esta regresión el modelo explica la totalidad de los datos utilizados para la construcción de la curva. A su vez, este resultado indica que el 98& de la variabilidad entre el modelo y los datos experimentales está dada por la variabilidad debida a la regresión (Ferreira 2007).
Prueba de significancia de la regresión
La prueba de significancia de la regresión hace una prueba con las siguientes hipótesis:
\(H_0:\) La Hipótesis nula que hay que comprobar implica que todos los coeficientes de la regresión son iguales a cero. \(H_1:\) La Hipótesis alternativa que se acepta cuando la nula es falsa, implica que al menos uno de los coeficientes es significativamente diferente de cero.
Esta prueba se hace al calcular el valor de F, y se rechaza la hipótesis nula a un nivel de signigicancia \(\alpha = 0,05\). con un valor de \(p = 1.104e-07\) confirma la significancia de almenos uno de los coeficientes, en este caso la pendiente.
Ajuste sin los niveles 75 mg Si /L y 100 mg Si / L de la curva
Construcción del modelo
Debido a que se obtuvieron dos resultados de absorbancia mayores a 1:
curva[c(9,8),]## Si.mg.L Absorbancia.u.a.
## 9 100 1.215
## 8 75 1.107Podemos revisar el comportamiento de la regresión restando estos dos puntos.
curva2 <- curva[-c(8,9),]
modelo2 <- lm( Absorbancia.u.a. ~ Si.mg.L , data = curva2)
summary(modelo2)##
## Call:
## lm(formula = Absorbancia.u.a. ~ Si.mg.L, data = curva2)
##
## Residuals:
## 1 2 3 4 5 6 7
## -0.0133425 -0.0127471 -0.0121516 -0.0003654 0.0416118 0.0075660 -0.0105712
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.0133425 0.0103684 1.287 0.255
## Si.mg.L 0.0124046 0.0004984 24.891 1.95e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.02186 on 5 degrees of freedom
## Multiple R-squared: 0.992, Adjusted R-squared: 0.9904
## F-statistic: 619.6 on 1 and 5 DF, p-value: 1.953e-06Observando la nueva línea de regresión y los datos expermientales:
library(ggplot2)
g <- ggplot(curva2, aes(Si.mg.L, Absorbancia.u.a.))
g + geom_point() + geom_smooth(method ="lm")## `geom_smooth()` using formula 'y ~ x'
Análisis de residuales para confirmar supuesto de homocedasticidad
Podemos realizar las mismas pruebas que para el primer modelo:
residuales2 <- modelo2$residuals
boxplot(residuales2)
abline(h = median(residuales2), col = "red", lty = 2)
median(residuales2)## [1] -0.01057122En este caso observamos que hay solo un dato atípico y que es posible gracias a la distribución dibujar los bigotes.
library(ggplot2)
g <- qplot(residuales2)
g + geom_density(alpha = 0.4,)## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
En este caso el grafico de densidad nos presenta una oblicuidad hacia la derecha, la cual podemos confirmar calculando el momento central:
library(moments)
kurtosis(residuales2)## [1] 3.79515skewness(residuales2)## [1] 1.474363En donde un valor de curtosis cercano a 3 indica que estamos más cerca de una distribución normal, con una oblicuidad positiva (cola hacia la derecha).
Podemos revisar la variación de los residuales en fución de las predicciones de absorbancia:
y.pred <- (curva2$Si.mg.L*modelo$coefficients[2]) + modelo$coefficients[1]
residuales <- data.frame(residuales = modelo2$residuals, y.pred = y.pred)
ggplot(residuales, aes(y.pred, residuales)) + geom_point()
Esta variación debería ser constante a lo largo de la variable predicha (absorbancia). en este gráfico de dispersión la variación aumenta a medida que aumenta la absorbancia, lo que indica que el error en x es significativo, contrario al supuesto de homocedasticidad.
Pruebas de parámetros de regresión
Para reportar los coeficientes de la regresión con intervalos de confianza:
Q <- qnorm(1-0.05/2)
intervalo_intercepto2 <- 0.0103684*Q
intervalo_pendiente2 <-0.0004984*Q
print(round(c(modelo2$coefficients[1],modelo2$coefficients[2]),3))## (Intercept) Si.mg.L
## 0.013 0.012print(round(c(intervalo_intercepto2,intervalo_pendiente2),3))## [1] 0.020 0.001Teniendo de nuevo el resultado de que el intercepto \(b_0 = 0.013 \pm 0.020 \; (p=0.255 )\) no es sigificativo para la regresión, y la pendiente \(b = 0.012\pm 0.001 \; (p=1.95\cdot10^{-6} )\) si es significativa.
De manera análoga al primer modelo se pueden evaluar los demás parámetros y pruebas de esta regresión.
Cuantificación de muestras problema
Se realizó la cuantificación de silicio en los extractos de digestión alcalina asistida por microondas y en la dilución preparada de silicio soluble en las muestras AFM402 y AFM 12863 siendo estas tierras de diatomeas y un fertilizante de silicio soluble edapot respectivamente.
En esta primera aproximación, hubo un error en la preparación de la
curva de calibración, y la lectura de la absorbancia de las muestras fue
realizada el 01/03/2021, mientras que la curva de calibración fue medida
el 03/03/2021, la cual fue utilizada para construir los objetos
modelo y modelo2.
Se guardan los nombres de las muestras y las absorbancias:
nombres <- c("blanco procedimiento 1",
"MRC Dil 2-1",
"402-D50-1",
"402-D50-2",
"402-D50-3",
"blanco procedimiento 2",
"402-D50-3",
"MRC Dil 2-2",
"402-D50-4",
"blanco procedimiento 3",
"402-D50-5",
"402-D50-6",
"402-D50-7",
"402-D50-8",
"Si Soluble D10",
"Si Soluble D20",
"MRC-1",
"MRC-2")
absorbancia <- c(0.048,0.237,0.124,0.120,0.062,0.048,0.106,0.271,0.159,0.046,0.117,0.179,0.134,0.132,0.847,0.472,0.432,0.447)Se juntan los datos en un solo data.frame :
muestras <- data.frame(muestra = nombres, absorbancia= absorbancia)
muestras## muestra absorbancia
## 1 blanco procedimiento 1 0.048
## 2 MRC Dil 2-1 0.237
## 3 402-D50-1 0.124
## 4 402-D50-2 0.120
## 5 402-D50-3 0.062
## 6 blanco procedimiento 2 0.048
## 7 402-D50-3 0.106
## 8 MRC Dil 2-2 0.271
## 9 402-D50-4 0.159
## 10 blanco procedimiento 3 0.046
## 11 402-D50-5 0.117
## 12 402-D50-6 0.179
## 13 402-D50-7 0.134
## 14 402-D50-8 0.132
## 15 Si Soluble D10 0.847
## 16 Si Soluble D20 0.472
## 17 MRC-1 0.432
## 18 MRC-2 0.447Luego, utilizando los coeficientes de regresión hallamos las concentraciones en mg Si/L
muestras$concentracion.sin.intercepto <- round((muestras$absorbancia /modelo2$coefficients[2]),3)
muestras $concentracion.con.intercepto <- round((muestras$absorbancia -modelo2$coefficients[1])/modelo2$coefficients[2],3)
knitr::kable(muestras)| muestra | absorbancia | concentracion.sin.intercepto | concentracion.con.intercepto |
|---|---|---|---|
| blanco procedimiento 1 | 0.048 | 3.870 | 2.794 |
| MRC Dil 2-1 | 0.237 | 19.106 | 18.030 |
| 402-D50-1 | 0.124 | 9.996 | 8.921 |
| 402-D50-2 | 0.120 | 9.674 | 8.598 |
| 402-D50-3 | 0.062 | 4.998 | 3.923 |
| blanco procedimiento 2 | 0.048 | 3.870 | 2.794 |
| 402-D50-3 | 0.106 | 8.545 | 7.470 |
| MRC Dil 2-2 | 0.271 | 21.847 | 20.771 |
| 402-D50-4 | 0.159 | 12.818 | 11.742 |
| blanco procedimiento 3 | 0.046 | 3.708 | 2.633 |
| 402-D50-5 | 0.117 | 9.432 | 8.356 |
| 402-D50-6 | 0.179 | 14.430 | 13.355 |
| 402-D50-7 | 0.134 | 10.802 | 9.727 |
| 402-D50-8 | 0.132 | 10.641 | 9.566 |
| Si Soluble D10 | 0.847 | 68.281 | 67.206 |
| Si Soluble D20 | 0.472 | 38.050 | 36.975 |
| MRC-1 | 0.432 | 34.826 | 33.750 |
| MRC-2 | 0.447 | 36.035 | 34.959 |
Posteriormente podemos convertir estas concentraciones a % p/p de \(SiO_2\):
- Primero, los extractos directos:
FD.MCR1 <- round(1.803807262,5)
FD.MCR2 <- round(1.850335307,5)
masas.muestras <- c(0, 0.101/FD.MCR1, 0.1027, 0.1001, 0.1023, 0, 0.1054, 0.1034/FD.MCR2, 0.1031, 0, 0.1010, 0.1011, 0.1025, 0.1113, NA, NA,0.101,0.1034
)
muestras$porc.SiO2 <- round(muestras$concentracion.sin.intercepto *(0.025*60.8*100)/(1000*28.09*masas.muestras), 4)
knitr::kable(muestras)| muestra | absorbancia | concentracion.sin.intercepto | concentracion.con.intercepto | porc.SiO2 |
|---|---|---|---|---|
| blanco procedimiento 1 | 0.048 | 3.870 | 2.794 | Inf |
| MRC Dil 2-1 | 0.237 | 19.106 | 18.030 | 1.8464 |
| 402-D50-1 | 0.124 | 9.996 | 8.921 | 0.5267 |
| 402-D50-2 | 0.120 | 9.674 | 8.598 | 0.5230 |
| 402-D50-3 | 0.062 | 4.998 | 3.923 | 0.2644 |
| blanco procedimiento 2 | 0.048 | 3.870 | 2.794 | Inf |
| 402-D50-3 | 0.106 | 8.545 | 7.470 | 0.4387 |
| MRC Dil 2-2 | 0.271 | 21.847 | 20.771 | 2.1155 |
| 402-D50-4 | 0.159 | 12.818 | 11.742 | 0.6727 |
| blanco procedimiento 3 | 0.046 | 3.708 | 2.633 | Inf |
| 402-D50-5 | 0.117 | 9.432 | 8.356 | 0.5053 |
| 402-D50-6 | 0.179 | 14.430 | 13.355 | 0.7723 |
| 402-D50-7 | 0.134 | 10.802 | 9.727 | 0.5703 |
| 402-D50-8 | 0.132 | 10.641 | 9.566 | 0.5173 |
| Si Soluble D10 | 0.847 | 68.281 | 67.206 | NA |
| Si Soluble D20 | 0.472 | 38.050 | 36.975 | NA |
| MRC-1 | 0.432 | 34.826 | 33.750 | 1.8658 |
| MRC-2 | 0.447 | 36.035 | 34.959 | 1.8858 |
Luego, las diluciones volumetricas de la muestra de silicio soluble en agua:
Vol.alic.D10.Si.Sol <- 1
Vol.alic.D20.Si.Sol <- 0.5
muestras$porc.SiO2[15] <- round( muestras$concentracion.con.intercepto[15]*(0.01*250*1*60.8*100)/(Vol.alic.D10.Si.Sol*1000*28.09*1.0182) ,4)
muestras$porc.SiO2[16] <- round( muestras$concentracion.con.intercepto[15]*(0.01*250*1*60.8*100)/(Vol.alic.D20.Si.Sol*1000*28.09*1.0182) ,4)
knitr::kable(muestras)| muestra | absorbancia | concentracion.sin.intercepto | concentracion.con.intercepto | porc.SiO2 |
|---|---|---|---|---|
| blanco procedimiento 1 | 0.048 | 3.870 | 2.794 | Inf |
| MRC Dil 2-1 | 0.237 | 19.106 | 18.030 | 1.8464 |
| 402-D50-1 | 0.124 | 9.996 | 8.921 | 0.5267 |
| 402-D50-2 | 0.120 | 9.674 | 8.598 | 0.5230 |
| 402-D50-3 | 0.062 | 4.998 | 3.923 | 0.2644 |
| blanco procedimiento 2 | 0.048 | 3.870 | 2.794 | Inf |
| 402-D50-3 | 0.106 | 8.545 | 7.470 | 0.4387 |
| MRC Dil 2-2 | 0.271 | 21.847 | 20.771 | 2.1155 |
| 402-D50-4 | 0.159 | 12.818 | 11.742 | 0.6727 |
| blanco procedimiento 3 | 0.046 | 3.708 | 2.633 | Inf |
| 402-D50-5 | 0.117 | 9.432 | 8.356 | 0.5053 |
| 402-D50-6 | 0.179 | 14.430 | 13.355 | 0.7723 |
| 402-D50-7 | 0.134 | 10.802 | 9.727 | 0.5703 |
| 402-D50-8 | 0.132 | 10.641 | 9.566 | 0.5173 |
| Si Soluble D10 | 0.847 | 68.281 | 67.206 | 35.7163 |
| Si Soluble D20 | 0.472 | 38.050 | 36.975 | 71.4327 |
| MRC-1 | 0.432 | 34.826 | 33.750 | 1.8658 |
| MRC-2 | 0.447 | 36.035 | 34.959 | 1.8858 |
En la tabla anterior están los resultados en \(\% \; SiO_2\):
- Para la muestra Muestra 402 podemos evaluar la distribución de los resultados de las réplicas, y verificar si podemos reportar el promedio como medida de tendencia central:
m.402 <- muestras[c(3,4,5,7,9,11:14),5]
boxplot(m.402)
En el gráfico de boxplot para los resultados en \(%SiO_2\) se encuentran 3 datos atípicos, las réplicas 3, 4 y 6. al removerlas podemos evaluar de nuevo la distribución sin efectos de escala por datos atípicos:
m.402.sin.atipicos <- m.402[c(3,4,6)]
boxplot(m.402.sin.atipicos)
Una vez evaluamos la distribución de los resultados sin datos atípicos, desde el boxplot se puede inferir que la muestra tiene oblicuidad hacia la derecha, lo que hará que los valores de la media y la mediana difieran:
summary(m.402.sin.atipicos)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.2644 0.3515 0.4387 0.4028 0.4720 0.5053Podemos confirmar la afirmación de oblicuidad negativa calculando el momento central de oblicuidad:
library(moments)
skewness(m.402.sin.atipicos)## [1] -0.4860191El valor medio del porcentaje de \(SiO_2\) para la muestra de tierra de diatomeas \(0.402\) difiere del reportado en un 99.46% :
error.diatomeas <-
((74.2-0.402)/74.2)*100
error.diatomeas ## [1] 99.45822Esto debido a la alta concentración de silicio en la muestra, la cual generó en el extracto una precipitación en forma de gel, y que a no logró pasar a través del filtro en el proceso de extracción. Este resultado permite proponer evaluar el rango de concentraciones de silicio en las muestras a trabajar, en el cual se mantenga un sesgo debido a la pérdida por precipitación en lo mínimo posible. Lo anterior podría lograrse disminuyendo la concentración de muestra en el extracto, o aumentando los reactivos. Otras variables podrían evaluarse para optimizar el proceso de extracción mediante un diseño experimental, similar al utilizado por Vereda Alonso 2011 variando el tiempo de digestión, la temperatura de digestión y las proporciones \(NaOH:muestra\) y \(H_2O_2:muestra\).
Para la muestra de silicio soluble se espera una concentración de ~ \(28 \%\)
tabla.Si.Soluble <- data.frame(muestra = muestras$muestra[c(15,16)], SiO2 = muestras[c(15,16),5])
tabla.Si.Soluble## muestra SiO2
## 1 Si Soluble D10 35.7163
## 2 Si Soluble D20 71.4327Valores que incluso en dilución, son extremadamente altos comparado a lo esperado, Lo cual puede deberse a la toma de una muestra no representativa del producto, y una dilución mayor a la adecuada en el momento de preparar el extracto.
Para el material de referencia certificado de roca fosfórica marroquí, el cual reporta un porcentaje en masa de \(SiO_2\) de \(2.01\%\).
tabla.MRC <- data.frame(muestra = muestras$muestra[c(2,8,17,18)],
SiO2 = muestras[c(2,8,17,18),5])
tabla.MRC$Error <- round(abs(((tabla.MRC$SiO2-2.01)/2.01)*100),3)
tabla.MRC## muestra SiO2 Error
## 1 MRC Dil 2-1 1.8464 8.139
## 2 MRC Dil 2-2 2.1155 5.249
## 3 MRC-1 1.8658 7.174
## 4 MRC-2 1.8858 6.179En la tabla anterior se presentan los resultados de la determinación del contenido de silicio en el material de referencia certificado de roca fosfórica, junto con el porcentaje de recuperación para cada réplica, en los cuales el máximo porcentaje de error relativo es \(8.14\%\). El comportamiento de los resultados para el material de referencia, soporta las afirmaciones sobre los resultados de las muestras de tierra de diatomeas y de silicio soluble, sobre un mayor efecto del proceso de extracción en el sesgo los resultados, y no del método de cuantificación.
Experimento 05/03/2021 Extractos paz de río - comité de silicio
En este experimento se analizaron un conjunto de extractos con ácido acético 0,514 M nombradas de la forma: B (Blanco), E (Extracto), EA (Extracto con adición de silicio) ,BA (Blanco con adición de silicio), las cuales se analizaron por duplicado para tener en cuenta el error de las diluciones volumétricas.
blanco.curva <- 0.01
concentraciones <- c(1,2,5,10,20,50,100)
absorbancias <- c(mean(0.016,0.013,0.016),
mean(0.021,0.020,0.020),
mean(0.039,0.041,0.043),
mean(0.083,0.085,0.085),
mean(0.160,0.0160,0.0156),
mean(0.391,0.392,0.390),
mean(0.694,0.695,0.696)
)
absorbancias <- absorbancias - blanco.curva
curva3 <- data.frame(Si.mg.L = concentraciones, absorbancia.u.a.=absorbancias)modelo.05.03.2021 <- lm(absorbancia.u.a.~Si.mg.L , data = curva3)
summary(modelo.05.03.2021)##
## Call:
## lm(formula = absorbancia.u.a. ~ Si.mg.L, data = curva3)
##
## Residuals:
## 1 2 3 4 5 6 7
## -0.0046744 -0.0066317 -0.0095038 -0.0002905 0.0071361 0.0294157 -0.0154515
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.0037170 0.0078207 0.475 0.655
## Si.mg.L 0.0069573 0.0001813 38.381 2.26e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.01619 on 5 degrees of freedom
## Multiple R-squared: 0.9966, Adjusted R-squared: 0.9959
## F-statistic: 1473 on 1 and 5 DF, p-value: 2.262e-07Q <- qnorm(1-0.05/2)
intervalo_intercepto3 <- Q* 0.0078207
intervalo_pendiente3 <- Q* 0.0001813
print(round(c(0.0037170 ,0.0069573 ),3))## [1] 0.004 0.007print(round(c(intervalo_intercepto3,intervalo_pendiente3),4))## [1] 0.0153 0.0004Para este modelo, el intercepto, \(b_0 = 0.004 \pm 0.0153 \;(p = 0.655 )\) no es significativamente diferente de cero, mientras que la pendiente \(b = 0.007 \pm 0.0004 \; (p = 2.26 \cdot 10^{-7})\) sí es significativamente diferente de cero.
library(ggplot2)
g <- ggplot(curva3, aes(Si.mg.L,absorbancia.u.a.))
g + geom_smooth(method = "lm") + geom_point()## `geom_smooth()` using formula 'y ~ x'
residuales3 <- modelo.05.03.2021$residuals
boxplot(residuales3)
abline(h = median(residuales3), col = "red", lty = 2)
median(residuales3)## [1] -0.004674382library(ggplot2)
g <- qplot(residuales3)
g + geom_density(alpha = 0.4,)## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
library(moments)
kurtosis(residuales3)## [1] 3.334163skewness(residuales3)## [1] 1.163268muestras2 <- c("B7",
"B7-2",
"AFM-09722",
"AFM-09722-2",
"B.A.7",
"B.A.7-2",
"E7",
"E7-2",
"E.A.7",
"E.A.7-2")
Abs.muestras <- c(mean(0.032,0.030,0.026),
mean(0.030,0.030,0.027),
mean(0.079,0.085,0.073),
mean(0.080,0.083,0.078),
mean(0.133,0.127,0.131),
mean(0.131,0.124,0.124),
mean(0.267,0.275,0.278),
mean(0.288,0.273,0.285),
mean(0.348,0.342,0.432),
mean(0.353,0.354,0.350))
muestras2.1 <- data.frame(Muestra = muestras2, Absorbanciau.a. = Abs.muestras)
muestras2.1$C.sin.intercepto <- muestras2.1$Absorbanciau.a./modelo.05.03.2021$coefficients[2]
muestras2.1$C.con.intercepto <- (muestras2.1$Absorbanciau.a.- modelo.05.03.2021$coefficients[1])/modelo.05.03.2021$coefficients[2]
print(muestras2.1)## Muestra Absorbanciau.a. C.sin.intercepto C.con.intercepto
## 1 B7 0.032 4.599456 4.065195
## 2 B7-2 0.030 4.311990 3.777729
## 3 AFM-09722 0.079 11.354907 10.820646
## 4 AFM-09722-2 0.080 11.498640 10.964379
## 5 B.A.7 0.133 19.116489 18.582229
## 6 B.A.7-2 0.131 18.829023 18.294763
## 7 E7 0.267 38.376712 37.842451
## 8 E7-2 0.288 41.395105 40.860844
## 9 E.A.7 0.348 50.019085 49.484824
## 10 E.A.7-2 0.353 50.737750 50.203489La tabla anterior presenta una concentración de Silicio en el blanco de 4.07 mg Si/L para la primera replica y 3.78 mg Si/L para la segunda réplica. Esto se atribuye al contacto de los extractos con material de vidrio durante el proceso de extracción, efecto que se minimizará en el siguiente experimento sin material de vidrio.
Cambio en el proceso de extracción “08/03/2021”
El día 08 de marzo de 2021 se preparó la colorimetría para el conjunto de muestras analizado en la sección anterior cambiando el recipiente de vidrio en la etapa de extracción por recipíentes de plástico para evitar la señal del blanco de procedimiento.
Se prepararon soluciones a partir de la curva de calibración para detectar la señal asociada a silicio en \(\lambda_{máx}= 660\; nm\) mediante el espectrofotómetro pharo 300.
ConcentracionesC.4 <- c(0,1,2,5,10,20,50,100)
Abs.H2O.HPLC <- c(0.012,0.09,0.013)
AbsorbanciasC.4 <- c(mean(0.05,-0.06,0.04),
mean(0.006,0.011,0.008),
mean(0.023,0.022,0.024),
mean(0.052,0.056,0.057),
mean(0.097,0.099,0.098),
mean(0.169,0.167,0.169),
mean(0.364,0.369,0.366),
mean(0.703,0.701,0.705))
curva4 <- data.frame(Si.mg.L = ConcentracionesC.4, Absorbancia.u.a. = AbsorbanciasC.4)modelo4 <- lm(Absorbancia.u.a. ~ Si.mg.L, data = curva4)
intercepto4 <- modelo4$coefficients[1]
pendiente4 <- modelo4$coefficients[2]
summary(modelo4)##
## Call:
## lm(formula = Absorbancia.u.a. ~ Si.mg.L, data = curva4)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.0234869 -0.0069110 -0.0008745 0.0070509 0.0273359
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.022664 0.007185 3.155 0.0197 *
## Si.mg.L 0.006823 0.000178 38.325 2.11e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.01652 on 6 degrees of freedom
## Multiple R-squared: 0.9959, Adjusted R-squared: 0.9953
## F-statistic: 1469 on 1 and 6 DF, p-value: 2.107e-08boxplot(modelo4$residuals)
y.pred4 <- (curva4$Si.mg.L*modelo4$coefficients[2]) + modelo4$coefficients[1]
residuales4 <- data.frame(residuales = modelo4$residuals, y.pred = y.pred4)
ggplot(residuales4, aes(y.pred, residuales)) + geom_point()
intervalo_intercepto4 <- 0.007185*Q
intervalo_pendiente4 <- 0.000178*Q
print(round(c(0.022664 ,0.006823),3))## [1] 0.023 0.007print(round(c(intervalo_intercepto4,intervalo_pendiente4),4))## [1] 0.0141 0.0003library(ggplot2)
g <- ggplot(curva4, aes(Si.mg.L,Absorbancia.u.a.))
g + geom_smooth(method = 'lm') + geom_point()## `geom_smooth()` using formula 'y ~ x'
Mediante el análisis de regresión se obtienen los parámetros \(b_0 = 0.023 \pm0.0141 (p= 0.0197)\) y \(b = 0.007 \pm 0.0003 (p= 2.11e-08)\)
En paralelo se prepararon las soluciones correspondientes a las muestras suministradas por el área de fertilizantes junto a tres controles analíticos, y una nueva preparación de la muestra de silicio soluble edapot AFM 12863, la cual fue preparada pesando 2.0144g y llevando a un volumen de 100 mL. A partir de el extracto se preparó una dilución 100, tomando 1 mL de extracto y llevandolo a un volúmen final de 100 mL. La dilución final fue seleccionada para la determinación de Si en solución mediante colorimetría:
muestras4 <- data.frame(Muestra = c("Control 1.mg.L",
"Control 20.mg.L",
"Control 100.mg.L",
"Blanco.A.A.1",
"E.A.1",
"E.A.1 + Silicio",
"E.Est Silicio 1",
"AFM 9722 - 1",
"AFM 9722 - 2",
"AFM 12863 D100"
)
)
muestras4$Absorbancia.u.a. <- c(mean(0.016,0.018,0.020),
mean(0.163,0.153,0.165),
mean(0.707,0.696,0.703),
mean(0.015,0.016,0.014),
mean(0.211,0.200,0.199),
mean(0.219,0.220,0.230),
mean(0.143,0.141,0.146),
mean(0.010,0.017,0.016),
mean(0.008,0.024,0.025),
mean(0.213,0.212,0.209))
muestras4$Si.mg.L.Con.Intercepto <- (muestras4$Absorbancia.u.a - intercepto4) / pendiente4
muestras4$Si.mg.L.Sin.Intercepto <- (muestras4$Absorbancia.u.a) / pendiente4
muestras4$Si.mg.L.Con.Intercepto <- round(muestras4$Si.mg.L.Con.Intercepto,2)
muestras4$Si.mg.L.Sin.Intercepto <- round(muestras4$Si.mg.L.Sin.Intercepto,2)
muestras4$SiO2 <- rep(NA, 10)
muestras4$SiO2[10] <-round((27.90*60.08*100*0.1*100)/(1000*28.09*2.0144),2)knitr::kable(muestras4)| Muestra | Absorbancia.u.a. | Si.mg.L.Con.Intercepto | Si.mg.L.Sin.Intercepto | SiO2 |
|---|---|---|---|---|
| Control 1.mg.L | 0.016 | -0.98 | 2.35 | NA |
| Control 20.mg.L | 0.163 | 20.57 | 23.89 | NA |
| Control 100.mg.L | 0.707 | 100.30 | 103.62 | NA |
| Blanco.A.A.1 | 0.015 | -1.12 | 2.20 | NA |
| E.A.1 | 0.211 | 27.60 | 30.93 | NA |
| E.A.1 + Silicio | 0.219 | 28.78 | 32.10 | NA |
| E.Est Silicio 1 | 0.143 | 17.64 | 20.96 | NA |
| AFM 9722 - 1 | 0.010 | -1.86 | 1.47 | NA |
| AFM 9722 - 2 | 0.008 | -2.15 | 1.17 | NA |
| AFM 12863 D100 | 0.213 | 27.90 | 31.22 | 29.62 |
Mediante los resultados podemos evaluar la exactitud del método mediante el porcentaje de recuperación de las adiciones realizadas, la adición del blanco fue de \(16,02 \;mg/L\) y la de el extracto paz de rio de \(15,80\;mg/L\)
muestras4$recuperacion <- rep(NA,10)
muestras4$recuperacion[6] <- ((muestras4$Si.mg.L.Con.Intercepto[6]- muestras4$Si.mg.L.Con.Intercepto[5])/15.80)*100
muestras4$recuperacion[7] <- ((muestras4$Si.mg.L.Con.Intercepto[7]- muestras4$Si.mg.L.Con.Intercepto[4])/16.02)*100Podemos también calcular la recuperación para la determinación de las adiciones mediante absorción atómica:
Muestra.PZR.AA <- 25.90
blanco.AA <- -5.58
Muestra.PZR.Ad <- 30.65
blanco.AA.Ad <- 16.64recuperaciones <- data.frame(muestra = c('Blanco + Si','muestra Paz de Rio + Si'),
recuperacion.A.A. = c( (blanco.AA.Ad - blanco.AA)/15.80*100 , (Muestra.PZR.Ad - Muestra.PZR.AA)/15.80*100 ),
recuperacion.UV = round(c(muestras4[7,5],muestras4[6,5]),3)
)
knitr::kable(recuperaciones)| muestra | recuperacion.A.A. | recuperacion.UV |
|---|---|---|
| Blanco + Si | 140.63291 | NA |
| muestra Paz de Rio + Si | 30.06329 | NA |
En la tabla anterior se presentan los porcentajes de recuperación de las adiciones de silicio, resultados que permiten evidenciar que la señal es mayor que la esperada en los blanco mediante anmbas técnicas espectrofotométricas, y a su vez la muestra de paz de río presenta una señal menor a la esperada en ambas técnicas espectrofotométricas. Lo cual indica un posible efecto matriz de la muestra paz de río, lo cual puede evidenciarse en la técnica colorimétrica cuando hay presencia de fósforo y hierro(Govett 1961).
Respecto a la muestra de silicio soluble Edapot AFM 12863 podemos obtener el porcentaje expresado como \(\% SiO_2\):
muestras4$SiO2[10]## [1] 29.62Podemos observar las concentraciones y absorbancias de las muestras en la curva de calibración:
library(ggplot2)
g <- ggplot(curva4, aes(Si.mg.L,Absorbancia.u.a.))
g + geom_smooth(method = 'lm') + geom_point() + geom_point(aes(Si.mg.L.Con.Intercepto,Absorbancia.u.a.), data = muestras4, colour = 'red')## `geom_smooth()` using formula 'y ~ x'
Debido a la alta dispersión en las respuestas de la curva de calibración a bajas concentraciones, existen predicciones para concentraciones negativas.
Conclusiones
Debido a la ausencia de una varianza constante en los residuales de las curvas de calibración, es aconsejable realizar diluciones gravimétricas para evitar el error en x.
Las mediciones a bajas concentraciones en el espectrofotómetro pharo 300 presentan una alta dispersión, por lo que es posible proponer una comparación de mediciones entre espectrofotómetros con el hach.
Es aconsejable estudiar el efecto matriz con la combinación de calibración externa, adición de estándar, y calibración de youden (Castells 2000).
Referencias
Saihua, L., Yunhe, X., Ji, X., Juan, H., Bocharnikova, E. A., & Matichenkov, V. V. (2018). Microwave digestion for colorimetric determination of total Si in plant and mineral samples. Communications in Soil Science and Plant Analysis, 49(7), 840-847.
Henderson, A. R. (2006). Testing experimental data for univariate normality. Clinica chimica acta, 366(1-2), 112-129.
Miller, J. N. (1991). Basic statistical methods for analytical chemistry. Part 2. Calibration and regression methods. A review. Analyst, 116(1), 3-14.
Mayerhöfer, T. G., Pahlow, S., & Popp, J. (2020). The Bouguer‐Beer‐Lambert law: Shining light on the obscure. ChemPhysChem. doi:10.1002/cphc.202000464
Alonso, E. V., de Torres, A. G., Cordero, M. S., & Pavón, J. C. (2011). Multivariate optimization of the synthesis and of the microwave dissolution of biomorphic silicon carbide ceramics. Microchemical Journal, 97(2), 101-108.
Ferreira, S. L. C., Bruns, R. E., da Silva, E. G. P., Dos Santos, W. N. L., Quintella, C. M., David, J. M., … & Neto, B. B. (2007). Statistical designs and response surface techniques for the optimization of chromatographic systems. Journal of chromatography A, 1158(1-2), 2-14.
Castells, R. C., & Castillo, M. A. (2000). Systematic errors: detection and correction by means of standard calibration, Youden calibration and standard additions method in conjunction with a method response model. Analytica chimica acta, 423(2), 179-185.
• Govett, G. J. S. (1961). Critical factors in the colorimetric determination of silica. Analytica Chimica Acta, 25(1), 69-80.