#Cargando librerias a utilizar

library(readxl)
library(PerformanceAnalytics)
library(ggplot2)
library(grid)
library(dplyr)
library(RColorBrewer)
library(scales)
library(gridExtra) 
library(ggthemes)
library(ggpubr)
#library(Rcmdr)
#library(RcmdrMisc)
library(NADA)
library(factoextra)
library(scales)
library(ggpubr)
library(tidyr)
library(cluster)
library(chemCal)
library(asbio)
library(ggridges)
library(tinytex)
library(knitr)
library(performance)
library(knitr)

Determinación de parámetros de desempeño relacionados con el modelo de regresión: linealidad,cedasticidad, Ld, Lc e incertidumbre

Esta regresión lineal se genera a partir de la repetición por duplicado de lectura de absorbancias de disoluciones con distintas concentraciones de dos analistas diferentes. Las curvas de calibracion con patrones de 0,1 a 0,5 mg/L permiten evaluar los datos correspondientes.

Determinación de linealidad

#Preparación de los datos

getwd()
## [1] "C:/Users/monge/OneDrive/Escritorio/Lab de Métodos de Análisis/Curso Métodos/Tarea Nitrito"
nitlinealidad<- read_excel("datavalidacion.xlsx",sheet = "linealidad")
View(nitlinealidad)
str(nitlinealidad)
## tibble [24 × 3] (S3: tbl_df/tbl/data.frame)
##  $ Patron       : chr [1:24] "Blanco" "1" "2" "3" ...
##  $ Concentracion: num [1:24] 0 0.101 0.202 0.303 0.404 ...
##  $ Absorbancia  : num [1:24] 0.001 0.079 0.166 0.261 0.358 0.449 0.001 0.079 0.169 0.263 ...
nitlinealidad$Concentracion<-as.numeric(nitlinealidad$Concentracion)
nitlinealidad$Absorbancia<-as.numeric(nitlinealidad$Absorbancia)
nitlinealidad$Patron<-as.factor(nitlinealidad$Patron)
regresion<-lm(Absorbancia~Concentracion,nitlinealidad) 

summary(regresion) 
## 
## Call:
## lm(formula = Absorbancia ~ Concentracion, data = nitlinealidad)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -0.0113 -0.0074 -0.0001  0.0045  0.0297 
## 
## Coefficients:
##                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   -0.004500   0.003164  -1.422    0.169    
## Concentracion  0.899359   0.010339  86.984   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.008743 on 22 degrees of freedom
## Multiple R-squared:  0.9971, Adjusted R-squared:  0.997 
## F-statistic:  7566 on 1 and 22 DF,  p-value: < 2.2e-16
ggplotRegression <- function (regresion) 
  
  require(ggplot2)
  
  ggplot(regresion$model, aes_string(x = names(regresion$model)[2], y = names(regresion$model)[1])) + 
    geom_point() +
    stat_smooth(method = "lm", col = "blue") +
    scale_y_continuous(labels = scales::number_format(accuracy = 0.001,decimal.mark = '.'))+
    scale_x_continuous(labels = scales::number_format(accuracy = 0.0001,decimal.mark = '.'))+
    labs(title = "Regresión lineal: validación del método espectrofotométrico",
         x= "Concentración (mg/L)", y="Absorbancia (UA)",
         subtitle = paste(" Pendiente =",signif(regresion$coef[[2]], 4),
                       "Intercepto =",signif(regresion$coef[[1]],4 ),
                       "r2=",signif(summary(regresion)$adj.r.squared, 4),
                       " p valor =",signif(summary(regresion)$coef[2,4], 3)))

plotlinealidad<-ggplotRegression(lm(Absorbancia~Concentracion, data = nitlinealidad)) 
plotlinealidad
## [1] TRUE

Análisis de residuales para comprobación de supuestos estadísticos

Aleatoriedad de residuales

plot(regresion$residuals) 

#gráfico de los residuales

Sumatoria cercana a cero

sum(regresion$residuals)# valor de la sumatoria de residuales
## [1] -1.290201e-17

Comprobación de la normalidad de los residuales mediante mediante prueba Shapiro-Wilk

shapiro.test(regresion$residuals) ### comprobacion de la normalidad de residuales
## 
##  Shapiro-Wilk normality test
## 
## data:  regresion$residuals
## W = 0.83359, p-value = 0.001101

Determinación de la cedasticidad

La evaluación de la cedasticidad en un enfoque espectrofotométrico implica analizar si la varianza de los datos es consistente en todas las concentraciones estudiadas. Por lo que se analiza cada patrón generado junto con sus varianzas para la determinación de la cedasticidad.

  ### convierte al patrón como un factor
nitlinealidad$Patron<-as.factor(nitlinealidad$Patron) 
### elimina los blancos para el análisis de cedasticidad
cedasticidad<-nitlinealidad %>% filter(Absorbancia > 0.01) 
View(cedasticidad)


bartlett.test(Absorbancia~Patron,cedasticidad) ### prueba Bartlett, 
## 
##  Bartlett test of homogeneity of variances
## 
## data:  Absorbancia by Patron
## Bartlett's K-squared = 21.623, df = 4, p-value = 0.0002382
plotcedasticidad<- ggplot(data = cedasticidad, aes(x = Patron, y = Absorbancia, colour = Patron)) +  ### grafico de cajas para observar la distribucion de los valores de absorbancia
  geom_jitter(width=0.20, size=3)+theme(legend.position = "none")+
  stat_summary(fun.data="mean_sdl", fun.args = list(mult=1), 
                 geom="errorbar", width = 0.2,color="black", 
                 position = position_dodge(0.8))+
  stat_summary(fun.data="mean_sdl", fun.args = list(mult=1), 
                 geom="pointrange", width = 0.6,color="black", 
                 position = position_dodge(0.8))+
  scale_y_continuous(labels = scales::number_format(accuracy = 0.001,decimal.mark = '.'))+
  labs(title = "Dispersión de los datos de absorbancia en los niveles calibración",
        x= "Nivel de concentración para la calibración", y="Absorbancia (UA)")
## Warning in stat_summary(fun.data = "mean_sdl", fun.args = list(mult = 1), :
## Ignoring unknown parameters: `width`
plotcedasticidad

Determinación de la incertidumbre por regresión lineal

La determinación de la incertidumbre se genera a partir de las muestras de absorbancia obtenidas de un MRC, y posteriormente calculando el intervalo de confianza del mismo.

muestra<-c(0.295,0.296,0.300, 0.309, 0.308, 0.313, 0.312) 
U<-inverse.predict(regresion,muestra,alpha = 0.05) ### determina la concentración, error
#estándar, IC y los límites de la predicción
U
## $Prediction
## [1] 0.3438164
## 
## $`Standard Error`
## [1] 0.004305511
## 
## $Confidence
## [1] 0.008929083
## 
## $`Confidence Limits`
## [1] 0.3348873 0.3527454
Ucal<-round(U$`Confidence`,digits=3)
Ucal #incertidumbre al 95% de confianza
## [1] 0.009

Determinación de los límites de detección y cuantificación

Se utiliza el método Miller & Miller el cual hace uso de la desviación estándar y el intercepto de la pendiente para calcular los límites de detección y cuantificación.

LoD<- round(((summary(regresion)$coef[[3]])*3.3)/(signif(regresion$coef[[2]], 4)), digits=3)
LoQ <- round(((summary(regresion)$coef[[3]])*10)/(signif(regresion$coef[[2]], 4)), digits=3)

LoD
## [1] 0.012
LoQ
## [1] 0.035

Determinación de repetibilidad y presición intermedia

Para analizar la repetibilidad y precisión intermedia de los resultados entre los dos estudiantes, se llevó a cabo una comparación entre los resultados de promedio, desviación estándar, los cuartiles (0, 0.25, 0.5, 0.75, 1), así como los valores mínimos y máximos.

Estadísitica descriptiva de los resultados

#Preparación de los datos
nitpres<- read_excel("datavalidacion.xlsx",sheet = "presicion")
View(nitpres)
 
nitpres$Estudiante<-as.factor(nitpres$Estudiante)
library(RcmdrMisc)
## Warning: package 'RcmdrMisc' was built under R version 4.3.3
## Loading required package: car
## Warning: package 'car' was built under R version 4.3.3
## Loading required package: carData
## Warning: package 'carData' was built under R version 4.3.3
## 
## Attaching package: 'car'
## The following object is masked from 'package:dplyr':
## 
##     recode
## Loading required package: sandwich
## Warning: package 'sandwich' was built under R version 4.3.3
descriptiva<-numSummary(nitpres[, c("Concentracion"), drop = FALSE], 
           groups = nitpres$Estudiante, statistics = c("mean", "sd", "IQR", "quantiles"),
           quantiles = c(0, 0.25, 0.5, 0.75, 1)) 

descriptiva
##        mean          sd    IQR    0%   25%   50%    75%  100% Concentracion:n
## D 0.3047143 0.007565586 0.0125 0.295 0.298 0.308 0.3105 0.313               7
## R 0.3162857 0.008769536 0.0130 0.305 0.310 0.315 0.3230 0.328               7

Repetibilidad para cada estudiante

Determinación del coeficiente de variación

CV_E1<-(descriptiva$table[1,2]/descriptiva$table[1,1])*100
CV_E1
## [1] 2.482846
CV_E2<-(descriptiva$table[2,2]/descriptiva$table[2,1])*100
CV_E2
## [1] 2.772663

Presición intermedia

Determinación de igualdad de varianzas por parte de ambos estudiantes

1.Análisis visual

plotpresi<-ggplot(data = nitpres, aes(x = Estudiante, y = Concentracion, colour = Estudiante)) +
  geom_boxplot(aes(),width=0.3) +geom_jitter(aes(shape=Estudiante),width=0.1,size=2) +
  scale_y_continuous(labels = scales::number_format(accuracy = 0.001,decimal.mark = '.'))+
  labs(title = "Dispersión de los datos de concentración",
       subtitle= "Factor de evaluación de la presición intermedia: Estudiante", y="Concentración (mg/L)")+
  theme(legend.position = "none") 
plotpresi ## grafico de cajas para visualizar la dispersiÓn de los datos de cada estudiante

2.Prueba de Homogeneidad de varianzas

La prueba F es muy potente para discriminar varianzas significativamente distintas, sin embargo posee la limitante que se supone que los datos deben de proveenir de una población normal.

Alternativamente, puede utilizarse la prueba de Levene la cuál presenta varias ventajas permite comparar 2 o más poblaciones,permite eligir entre diversos estadísticos (mediana, media o media truncada) y además no requiere que los datos sean normales como si lo es para el caso de la prueba F de Fisher.

est1<-nitpres %>% filter(Estudiante == "D")#filtrado de los datos por estudiante
est2<-nitpres %>% filter(Estudiante == "R")

shapiro.test(est1$Concentracion)### prueba de normalidad para los valores de concnetracion
## 
##  Shapiro-Wilk normality test
## 
## data:  est1$Concentracion
## W = 0.87495, p-value = 0.2049
shapiro.test(est2$Concentracion)
## 
##  Shapiro-Wilk normality test
## 
## data:  est2$Concentracion
## W = 0.94799, p-value = 0.7114
var.test(x = est1$Concentracion,y = est2$Concentracion)## prueba F para probar varianzas
## 
##  F test to compare two variances
## 
## data:  est1$Concentracion and est2$Concentracion
## F = 0.74427, num df = 6, denom df = 6, p-value = 0.729
## alternative hypothesis: true ratio of variances is not equal to 1
## 95 percent confidence interval:
##  0.1278872 4.3314845
## sample estimates:
## ratio of variances 
##          0.7442724
### 
leveneTest(y = nitpres$Concentracion, group = nitpres$Estudiante, center = "mean") 
## Levene's Test for Homogeneity of Variance (center = "mean")
##       Df F value Pr(>F)
## group  1  0.0997 0.7576
##       12

Según los resultados de cualesquiera de las dos pruebas, se observa que no existe diferencia estadísticamente significativa en las varianzas de ambos estudiantes, con esto sevalida la presición intermedia.

  1. Comprobación de igualdad de medias

El análisis puede complementarse con una prueba t student para analizar si hay diferencia significativa en el valor medio obtenido por ambos estudiantes

### Diferencia de medias
t.test(est1$Concentracion,est2$Concentracion,mu=0,paired=TRUE,var.equal=TRUE)
## 
##  Paired t-test
## 
## data:  est1$Concentracion and est2$Concentracion
## t = -5.2541, df = 6, p-value = 0.001912
## alternative hypothesis: true mean difference is not equal to 0
## 95 percent confidence interval:
##  -0.016960382 -0.006182475
## sample estimates:
## mean difference 
##     -0.01157143

Determinación de la Veracidad

Para la determinación de veracidad se contraste el valor medio obtenido contra el valor de referencia aportado por el MRC mediante el comando t.test

t.test(est1$Concentracion,mu=0.3505,alternative = "less")
## 
##  One Sample t-test
## 
## data:  est1$Concentracion
## t = -16.012, df = 6, p-value = 1.885e-06
## alternative hypothesis: true mean is less than 0.3505
## 95 percent confidence interval:
##       -Inf 0.3102709
## sample estimates:
## mean of x 
## 0.3047143
#Se realiza una prueba de una cola, el estudiante 1 obtuvo una media de 0.3047. 
#Hipótesis nula: la media del estudiantes 1 es igual a 0.3505 (valor del MRC) 
#hipótesis alterna: La media del estudiante 1 es menor a 0.3505

t.test(est2$Concentracion,mu=0.3505,alternative = "less") 
## 
##  One Sample t-test
## 
## data:  est2$Concentracion
## t = -10.322, df = 6, p-value = 2.416e-05
## alternative hypothesis: true mean is less than 0.3505
## 95 percent confidence interval:
##       -Inf 0.3227265
## sample estimates:
## mean of x 
## 0.3162857
#Se realiza una prueba de una cola, el estudiante 2 tuvo una media igual a 0.3163
#Hipótesis nula: la media del estudiantes 1 es igual a 0.3505 (valor del MRC)
#Hipótesis alterna: La media del estudiante 1 es MAYOR a 0.3505

## NOTA: probar que la media no difiere del MRC es solo una parte de la prueba de veracidad