if(!require(pacman)) install.packages("pacman", dependencies=TRUE)
## Loading required package: pacman
## Warning: package 'pacman' was built under R version 4.2.3
library("pacman")
p_load("vroom",
"ggplot2",
"ggrepel",
"dplyr",
"tidyverse")
Datos <- vroom(file = "https://raw.githubusercontent.com/ManuelLaraMVZ/resultados_PCR_practica/refs/heads/main/Cts.csv")
## Rows: 29 Columns: 5
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (4): Well, Grupo, Practica, Fluor
## dbl (1): Cq
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Datos
## # A tibble: 29 × 5
## Well Grupo Practica Fluor Cq
## <chr> <chr> <chr> <chr> <dbl>
## 1 A01 G1 Relativa SYBR 30.8
## 2 B01 G1 Relativa SYBR 41
## 3 C01 G1 Relativa SYBR 26.1
## 4 D01 G1 Relativa SYBR 41
## 5 E01 G1 Relativa SYBR 41
## 6 F01 G1 Relativa SYBR 22.3
## 7 A01 G2 Relativa SYBR 30.7
## 8 B01 G2 Relativa SYBR 20.3
## 9 C01 G2 Relativa SYBR 18.7
## 10 D01 G2 Relativa SYBR 20.5
## # … with 19 more rows
Seleccion de datos de la curva
datos_rt_pcr <- Datos %>%
filter(Grupo == "Curva") %>%
select(1,5) %>%
mutate(Well = as.numeric(Well),
Curva = Well,
logConc = log10(Curva)) %>%
select(3,4,2)
datos_rt_pcr
## # A tibble: 5 × 3
## Curva logConc Cq
## <dbl> <dbl> <dbl>
## 1 100 2 11.2
## 2 10 1 15.2
## 3 1 0 17.6
## 4 0.1 -1 20.1
## 5 0.01 -2 23.9
grafica <- ggplot(mapping = aes(x = datos_rt_pcr$logConc, y = datos_rt_pcr$Cq)) +
geom_point(color = "#900C3F", size = 5) +
theme_classic()
grafica
REalizar cálculos para el ajuste lineal
# Ajustar modelo de regresión lineal
modelo_pcr <- lm(Cq ~ logConc, data = datos_rt_pcr)
modelo_pcr
##
## Call:
## lm(formula = Cq ~ logConc, data = datos_rt_pcr)
##
## Coefficients:
## (Intercept) logConc
## 17.608 -3.024
# Obtener coeficientes de la ecuación
coeficientes <- coef(modelo_pcr)
pendiente <- round(coeficientes[2], 3)
intercepto <- round(coeficientes[1], 3)
coeficientes
## (Intercept) logConc
## 17.608002 -3.024013
pendiente
## logConc
## -3.024
intercepto
## (Intercept)
## 17.608
# Crear puntos de predicción para la línea ajustada
predicciones <- data.frame(logConc = seq(min(datos_rt_pcr$logConc), max(datos_rt_pcr$logConc), length.out = 100))
predicciones$Cq <- predict(modelo_pcr, newdata = predicciones)
predicciones
## logConc Cq
## 1 -2.00000000 23.65603
## 2 -1.95959596 23.53385
## 3 -1.91919192 23.41166
## 4 -1.87878788 23.28948
## 5 -1.83838384 23.16730
## 6 -1.79797980 23.04512
## 7 -1.75757576 22.92293
## 8 -1.71717172 22.80075
## 9 -1.67676768 22.67857
## 10 -1.63636364 22.55639
## 11 -1.59595960 22.43420
## 12 -1.55555556 22.31202
## 13 -1.51515152 22.18984
## 14 -1.47474747 22.06766
## 15 -1.43434343 21.94548
## 16 -1.39393939 21.82329
## 17 -1.35353535 21.70111
## 18 -1.31313131 21.57893
## 19 -1.27272727 21.45675
## 20 -1.23232323 21.33456
## 21 -1.19191919 21.21238
## 22 -1.15151515 21.09020
## 23 -1.11111111 20.96802
## 24 -1.07070707 20.84583
## 25 -1.03030303 20.72365
## 26 -0.98989899 20.60147
## 27 -0.94949495 20.47929
## 28 -0.90909091 20.35710
## 29 -0.86868687 20.23492
## 30 -0.82828283 20.11274
## 31 -0.78787879 19.99056
## 32 -0.74747475 19.86838
## 33 -0.70707071 19.74619
## 34 -0.66666667 19.62401
## 35 -0.62626263 19.50183
## 36 -0.58585859 19.37965
## 37 -0.54545455 19.25746
## 38 -0.50505051 19.13528
## 39 -0.46464646 19.01310
## 40 -0.42424242 18.89092
## 41 -0.38383838 18.76873
## 42 -0.34343434 18.64655
## 43 -0.30303030 18.52437
## 44 -0.26262626 18.40219
## 45 -0.22222222 18.28000
## 46 -0.18181818 18.15782
## 47 -0.14141414 18.03564
## 48 -0.10101010 17.91346
## 49 -0.06060606 17.79128
## 50 -0.02020202 17.66909
## 51 0.02020202 17.54691
## 52 0.06060606 17.42473
## 53 0.10101010 17.30255
## 54 0.14141414 17.18036
## 55 0.18181818 17.05818
## 56 0.22222222 16.93600
## 57 0.26262626 16.81382
## 58 0.30303030 16.69163
## 59 0.34343434 16.56945
## 60 0.38383838 16.44727
## 61 0.42424242 16.32509
## 62 0.46464646 16.20290
## 63 0.50505051 16.08072
## 64 0.54545455 15.95854
## 65 0.58585859 15.83636
## 66 0.62626263 15.71418
## 67 0.66666667 15.59199
## 68 0.70707071 15.46981
## 69 0.74747475 15.34763
## 70 0.78787879 15.22545
## 71 0.82828283 15.10326
## 72 0.86868687 14.98108
## 73 0.90909091 14.85890
## 74 0.94949495 14.73672
## 75 0.98989899 14.61453
## 76 1.03030303 14.49235
## 77 1.07070707 14.37017
## 78 1.11111111 14.24799
## 79 1.15151515 14.12580
## 80 1.19191919 14.00362
## 81 1.23232323 13.88144
## 82 1.27272727 13.75926
## 83 1.31313131 13.63708
## 84 1.35353535 13.51489
## 85 1.39393939 13.39271
## 86 1.43434343 13.27053
## 87 1.47474747 13.14835
## 88 1.51515152 13.02616
## 89 1.55555556 12.90398
## 90 1.59595960 12.78180
## 91 1.63636364 12.65962
## 92 1.67676768 12.53743
## 93 1.71717172 12.41525
## 94 1.75757576 12.29307
## 95 1.79797980 12.17089
## 96 1.83838384 12.04870
## 97 1.87878788 11.92652
## 98 1.91919192 11.80434
## 99 1.95959596 11.68216
## 100 2.00000000 11.55998
Obtener el valor de ajute (R^2)
# Obtener resumen del modelo
resumen_modelo <- summary(modelo_pcr)
resumen_modelo
##
## Call:
## lm(formula = Cq ~ logConc, data = datos_rt_pcr)
##
## Residuals:
## 1 2 3 4 5
## -0.328525 0.611422 -0.000825 -0.518515 0.236443
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 17.6080 0.2319 75.94 5.03e-06 ***
## logConc -3.0240 0.1640 -18.44 0.000348 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.5185 on 3 degrees of freedom
## Multiple R-squared: 0.9913, Adjusted R-squared: 0.9883
## F-statistic: 340.1 on 1 and 3 DF, p-value: 0.0003478
Extraer la ecuación de la recta y el valor de R^2
# Extraer el valor de R²
R2 <- resumen_modelo$r.squared
# Crear la ecuación en formato de texto
ecuacion_texto <- paste0("y = ", pendiente, "x + ", intercepto, "\nR² = ", round(R2, 4))
# Mostrar el valor de R²
cat("y = ", pendiente, "x + ", intercepto, "\nR² = ", round(R2, 4),"\n")
## y = -3.024 x + 17.608
## R² = 0.9913
Graficación
grafica2 <- grafica + # Puntos de datos reales
geom_smooth(data = predicciones, aes(x = logConc, y = Cq), color = "#FFC300", size = 1.5) + # Ajuste lineal
labs(
title = "Curva Estándar de RT-PCR en Tiempo Real",
subtitle = "Simulación de datos de RT-PCR",
caption = "Diseñó: Manuel Lara",
x = "Log10(Concentración de ADN [pg/µL])",
y = "Cq (Ct)") +
theme_classic(base_size = 15) +
scale_y_continuous(labels = scales::number_format(accuracy = 1),) + # Redondeo de 3 dígitos en Y
theme(
plot.title = element_text(hjust = 0.5, face = "bold"),
axis.title.x = element_text(face = "bold"),
axis.title.y = element_text(face = "bold")
) +
#geom_smooth(method = "lm", se = FALSE, color="#FFC300") + # Línea de regresión
annotate("text",
x = max(datos_rt_pcr$logConc) -1.5, # Posición hacia la derecha
y = min(datos_rt_pcr$Cq) + 12, # Ajustar para que no se sobreponga
label = ecuacion_texto,
size = 5, color = "#581845", fontface = "bold", hjust = 0)
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
grafica2
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'
Análisis de las muestras:
Datos
## # A tibble: 29 × 5
## Well Grupo Practica Fluor Cq
## <chr> <chr> <chr> <chr> <dbl>
## 1 A01 G1 Relativa SYBR 30.8
## 2 B01 G1 Relativa SYBR 41
## 3 C01 G1 Relativa SYBR 26.1
## 4 D01 G1 Relativa SYBR 41
## 5 E01 G1 Relativa SYBR 41
## 6 F01 G1 Relativa SYBR 22.3
## 7 A01 G2 Relativa SYBR 30.7
## 8 B01 G2 Relativa SYBR 20.3
## 9 C01 G2 Relativa SYBR 18.7
## 10 D01 G2 Relativa SYBR 20.5
## # … with 19 more rows
Analisis <- Datos %>%
filter(Practica == "Absoluta",
Well == "Profesor",
Grupo == "G1") %>%
select(1,2,5)
Analisis
## # A tibble: 1 × 3
## Well Grupo Cq
## <chr> <chr> <dbl>
## 1 Profesor G1 18.6
Colocar los Ct obtenidos del profesor
# Ct (Cq) observado
Ct_Control <- Analisis$Cq
Ct_Control
## [1] 18.57611
# Calcular log10(Concentración)
logConc_calculado_Cx <- (Ct_Control - intercepto) / pendiente
logConc_calculado_Cx
## (Intercept)
## -0.3201429
# Obtener concentración real
Conc_calculada_Cx <- 10^logConc_calculado_Cx
Conc_calculada_Cx
## (Intercept)
## 0.4784726
# Crear un data frame con el nuevo punto
nuevo_punto_Cx <- data.frame(logConc = logConc_calculado_Cx, Cq = Ct_Control)
nuevo_punto_Cx
## logConc Cq
## (Intercept) -0.3201429 18.57611
# Mostrar resultado
cat("La concentración estimada para la muestra del profesor es:", round(Conc_calculada_Cx, 4), "pg/µL\n")
## La concentración estimada para la muestra del profesor es: 0.4785 pg/µL
Graficamos el punto en la recta
interpolacion <- grafica2+
geom_point(data = nuevo_punto_Cx,
aes(x = logConc, y = Cq),
color = "black",
size = 6, shape = 20) + # Nuevo punto
geom_segment(aes(x = min(datos_rt_pcr$logConc),
xend = logConc_calculado_Cx,
y = Ct_Control,
yend = Ct_Control),
linetype = "dotted",
color = "#2f2d2c",
size = 1.5) + # Línea horizontal hasta el punto
geom_segment(aes(x = logConc_calculado_Cx, xend = logConc_calculado_Cx,
y = min(datos_rt_pcr$Cq), yend = Ct_Control),
linetype = "dotted",
color = "#2f2d2c",
size = 1.5)+
annotate("text",
x = logConc_calculado_Cx,
y = Ct_Control - 1,
label = paste0(" [Prof]= ",
round(Conc_calculada_Cx, 4), " pg/µL"),
color = "black", size = 5,
fontface = "bold",
hjust = -0.007, vjust = -1.7)
interpolacion
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'
Datos
## # A tibble: 29 × 5
## Well Grupo Practica Fluor Cq
## <chr> <chr> <chr> <chr> <dbl>
## 1 A01 G1 Relativa SYBR 30.8
## 2 B01 G1 Relativa SYBR 41
## 3 C01 G1 Relativa SYBR 26.1
## 4 D01 G1 Relativa SYBR 41
## 5 E01 G1 Relativa SYBR 41
## 6 F01 G1 Relativa SYBR 22.3
## 7 A01 G2 Relativa SYBR 30.7
## 8 B01 G2 Relativa SYBR 20.3
## 9 C01 G2 Relativa SYBR 18.7
## 10 D01 G2 Relativa SYBR 20.5
## # … with 19 more rows
Analisis2 <- Datos %>%
filter(Practica == "Absoluta",
Well == "Profesor",
Grupo == "G2") %>%
select(1,2,5)
Analisis2
## # A tibble: 1 × 3
## Well Grupo Cq
## <chr> <chr> <dbl>
## 1 Profesor G2 20.9
Colocar los Ct obtenidos del equipo
# Ct (Cq) observado
Ct_equipo<- Analisis2$Cq
Ct_equipo
## [1] 20.90737
# Calcular log10(Concentración)
logConc_calculado_equipo <- (Ct_equipo - intercepto) / pendiente
logConc_calculado_equipo
## (Intercept)
## -1.091061
# Obtener concentración real
Conc_calculada_equipo <- 10^logConc_calculado_equipo
Conc_calculada_equipo
## (Intercept)
## 0.08108467
# Crear un data frame con el nuevo punto
nuevo_punto_equipo <- data.frame(logConc = logConc_calculado_equipo, Cq = Ct_equipo)
nuevo_punto_equipo
## logConc Cq
## (Intercept) -1.091061 20.90737
# Mostrar resultado
cat("La concentración estimada para la muestra del equipo es:", round(Conc_calculada_equipo, 4), "pg/µL\n")
## La concentración estimada para la muestra del equipo es: 0.0811 pg/µL
Interpolación en la gráfica
interpolacion2 <- interpolacion+
geom_point(data = nuevo_punto_equipo,
aes(x = logConc, y = Cq),
color = "black",
size = 6, shape = 20) + # Nuevo punto en verde
geom_segment(aes(x = min(datos_rt_pcr$logConc),
xend = logConc_calculado_equipo,
y = Ct_equipo,
yend = Ct_equipo),
linetype = "dotted",
color = "#2f2d2c",
size = 1.5) + # Línea horizontal hasta el punto
geom_segment(aes(x = logConc_calculado_equipo, xend = logConc_calculado_equipo,
y = min(datos_rt_pcr$Cq), yend = Ct_equipo),
linetype = "dotted",
color = "#2f2d2c",
size = 1.5)+
annotate("text",
x = logConc_calculado_equipo,
y = Ct_equipo - 1,
label = paste0(" [Equipo]= ",
round(Conc_calculada_equipo, 4), " pg/µL"),
color = "black", size = 5,
fontface = "bold",
hjust = -0.007, vjust = -1.7)
interpolacion2
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'
Comparación con respecto al profesor
Comparacion <- data.frame(Condicion=c("Profesor", "Equipo"),
Concentracion= c(Conc_calculada_Cx, Conc_calculada_equipo))
Comparacion
## Condicion Concentracion
## 1 Profesor 0.47847260
## 2 Equipo 0.08108467
# Calcular el máximo
max_conc <- max(Comparacion$Concentracion, na.rm = TRUE)
# Calcular incremento, asegurando que no sea 0
break_by <- round(max_conc / 10, 1)
break_by <- ifelse(break_by == 0, 0.1, break_by)
# Crear la gráfica
Barras <- ggplot(Comparacion, aes(x = Condicion, y = Concentracion, fill = Condicion)) +
geom_bar(stat = "identity", color = "black", width = 0.9) +
labs(title = "Comparación de Concentraciones",
x = "Condición",
y = "Concentración (pg/µL)") +
theme_classic(base_size = 14) +
scale_fill_manual(values = c("Profesor" = "#f8c471", "Equipo" = "#2980b9")) +
scale_y_continuous(limits = c(0, max_conc+ 0.1),
breaks = seq(0, max_conc + 1, by = 0.1)) +
theme(
axis.title = element_text(face = "bold", size = 16),
axis.text = element_text(face = "bold", size = 14),
axis.line = element_line(size = 1.5, color = "black"),
legend.position = "none")+
geom_text(aes(label = paste0(round(Concentracion, 2), " pg/µL")),
vjust = -0.5, size = 5)
## Warning: The `size` argument of `element_line()` is deprecated as of ggplot2 3.4.0.
## ℹ Please use the `linewidth` argument instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
Barras
Porcentaje <- Comparacion %>%
mutate(Porcentaje = (Concentracion ) *100/ Concentracion[Condicion == "Profesor"]) %>%
select(1,3)
Porcentaje
## Condicion Porcentaje
## 1 Profesor 100.00000
## 2 Equipo 16.94656
# Calcular el máximo valor
max_porcentaje <- max(Porcentaje$Porcentaje, na.rm = TRUE)
# Calcular incremento, asegurando que no sea 0
break_by <- round(max_porcentaje / 10, 0)
break_by <- ifelse(break_by == 0, 1, break_by) # Asegura que by no sea 0
# Crear la gráfica
Barras2 <- ggplot(Porcentaje, aes(x = Condicion, y = Porcentaje, fill = Condicion)) +
geom_bar(stat = "identity", color = "black", width = 0.9) +
labs(
title = "Comparación de Concentraciones",
x = "Condición",
y = "Tasa de cambio (u. a.)"
) +
theme_classic(base_size = 14) +
scale_fill_manual(values = c("Profesor" = "#f8c471",
"Equipo" = "#2980b9")) +
scale_y_continuous(limits = c(0, max_porcentaje +20),
breaks = seq(0, max_porcentaje + 20, by = 10)) +
theme(axis.title = element_text(face = "bold", size = 16),
axis.text = element_text(face = "bold", size = 14),
axis.line = element_line(size = 1.5, color = "black"),
legend.position = "none") +
geom_text(aes(label = paste0(round(Porcentaje, 2), " %")), # <- Aquí está la etiqueta con unidad
vjust = -0.5, size = 5)
# Mostrar la gráfica
Barras2