MULTI-LEVEL ANALYSIS OF STRESS RESPONSES IN AMPHIBIAN OF EXTREME HABITATS
Alexandra Delgadillo M.
12/8/2020
- Pregunta de investigación:
- Resultados:
- Resultados morfológicos
- Analisis para evaluar la repetibilidad de las medidas tomadas en cada renacuajo:
- Análisis de componentes pricipales variables transformadas
#MULTI-LEVEL ANALYSIS OF STRESS RESPONSES IN AMPHIBIAN OF EXTREME HABITATS
Resultados Salida de Campo 2019 - Afloramientos Rocosos Reserva Natural Bojonawi
Pregunta de investigación:
El nivel de desecación de las charcas afecta la morfología y el desempeño locomotor de renacuajos de Leptodactylus lithonaetes
Resultados:
Estos resultados están divididos en dos secciones:
Los resultados morfológicos.
Los resultados en el desempeño locomotor.
Resultados morfológicos
Se utilizaron 7 posturas diferentes. Cada postura se dividió en 3 tratamientos (HW-High water treatment, LW-low water treatment y DD-Dry down treatment. para un total de 21 cajas. En cada caja se introdujeron 10 renacuajos. Al final del experimento, de cada caja se tomaron entre 8 y 10 renacuajos y a cada renacuajo se le tomaron tres veces las siguientes medidas morfológicas:
Altura de la cabeza
Altura de la cola
Longitud de la cola
Perímetro del músculo
Longitud de la cabeza
*Eliminé la longitud total por estar relacionada directamente con la longitud de la cola y la longitud de la cabeza (evitar redundancia)
knitr::include_graphics("dis_morf.png")
Analisis para evaluar la repetibilidad de las medidas tomadas en cada renacuajo:
La Base de datos
datos_morfo <- read.table("error_prueba_completo_corregido.txt",h=T)
head(datos_morfo)## Id_Random_Foto Tratamiento Replica Caja Individuo Medicion Altura_Cabeza
## 1 1 AE R2 16 4 1 2.787
## 2 1 AE R2 16 4 2 2.308
## 3 1 AE R2 16 4 3 2.364
## 4 2 AE R3 17 1 1 2.183
## 5 2 AE R3 17 1 2 2.164
## 6 2 AE R3 17 1 3 2.146
## Longitud_Total Altura_Cola Longitud_Cola Perim_Musculo1 Longitud_Cabeza
## 1 25.214 1.614 17.883 38.806 7.633
## 2 24.261 1.416 16.859 37.148 7.723
## 3 25.454 1.666 17.839 38.997 7.543
## 4 23.777 1.474 16.800 34.641 7.115
## 5 23.795 1.455 16.853 34.516 7.124
## 6 23.749 1.474 16.882 34.992 7.104- Para evaluar la repetibilidad entre las tres medidas tomadas para cada animal se realizaron modelos mixtos NULOS sin varible independiente. A partir de ese modelo se calculó el valor ICC (intraclass correlation coefficient) que permite determinar la “estabilidad en la medida” (Tomado de Finch et al., 2014: Multilevel modeling using R-p44).
ICC>0.90: Excelente
ICC entre 0,75 y 0,89: Bueno
ICC entre 0,50 y 0,74: Moderado
ICC<0.50: Malo
Se usó el paquete performance para ese cálculo de ICC.
Librerias
library(Matrix)
library(lme4)
library(lmerTest)#para calcular valor de p en los modelos
library(performance)
library(sjPlot) ###para hacer las tablas####Long_cabeza:
mixtolong_ca<-lmer(Longitud_Cabeza~1+(1|Id_Random_Foto),data=datos_morfo)
summary(mixtolong_ca)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Longitud_Cabeza ~ 1 + (1 | Id_Random_Foto)
## Data: datos_morfo
##
## REML criterion at convergence: 329.6
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.0533 -0.1848 -0.0115 0.1830 6.0879
##
## Random effects:
## Groups Name Variance Std.Dev.
## Id_Random_Foto (Intercept) 0.9498 0.9746
## Residual 0.0213 0.1459
## Number of obs: 525, groups: Id_Random_Foto, 175
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 6.13404 0.07394 173.99998 82.95 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1icc(mixtolong_ca, by_group = TRUE)## # ICC by Group
##
## Group | ICC
## ----------------------
## Id_Random_Foto | 0.978tab_model(mixtolong_ca)| Longitud_Cabeza | |||
|---|---|---|---|
| Predictors | Estimates | CI | p |
| (Intercept) | 6.13 | 5.99 – 6.28 | <0.001 |
| Random Effects | |||
| σ2 | 0.02 | ||
| τ00 Id_Random_Foto | 0.95 | ||
| ICC | 0.98 | ||
| N Id_Random_Foto | 175 | ||
| Observations | 525 | ||
| Marginal R2 / Conditional R2 | 0.000 / 0.978 | ||
####Altura_Cabeza:
mixtoalt_cabeza<-lmer(Altura_Cabeza~1+(1|Id_Random_Foto),data=datos_morfo)
summary(mixtoalt_cabeza)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Altura_Cabeza ~ 1 + (1 | Id_Random_Foto)
## Data: datos_morfo
##
## REML criterion at convergence: -469.7
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.6323 -0.3158 -0.0283 0.2970 4.3991
##
## Random effects:
## Groups Name Variance Std.Dev.
## Id_Random_Foto (Intercept) 0.087848 0.29639
## Residual 0.007021 0.08379
## Number of obs: 525, groups: Id_Random_Foto, 175
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 2.0141 0.0227 174.0000 88.72 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1icc(mixtoalt_cabeza, by_group = TRUE)## # ICC by Group
##
## Group | ICC
## ----------------------
## Id_Random_Foto | 0.926tab_model(mixtoalt_cabeza)| Altura_Cabeza | |||
|---|---|---|---|
| Predictors | Estimates | CI | p |
| (Intercept) | 2.01 | 1.97 – 2.06 | <0.001 |
| Random Effects | |||
| σ2 | 0.01 | ||
| τ00 Id_Random_Foto | 0.09 | ||
| ICC | 0.93 | ||
| N Id_Random_Foto | 175 | ||
| Observations | 525 | ||
| Marginal R2 / Conditional R2 | 0.000 / 0.926 | ||
####Longitud_Total
mixtolong_total<-lmer(Longitud_Total~1+(1|Id_Random_Foto),data=datos_morfo)
summary(mixtolong_total)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Longitud_Total ~ 1 + (1 | Id_Random_Foto)
## Data: datos_morfo
##
## REML criterion at convergence: 1518.8
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.2728 -0.1337 -0.0034 0.1411 6.1981
##
## Random effects:
## Groups Name Variance Std.Dev.
## Id_Random_Foto (Intercept) 9.1547 3.0257
## Residual 0.2064 0.4543
## Number of obs: 525, groups: Id_Random_Foto, 175
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 20.5020 0.2296 174.0000 89.3 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1icc(mixtolong_total, by_group = TRUE)## # ICC by Group
##
## Group | ICC
## ----------------------
## Id_Random_Foto | 0.978tab_model(mixtolong_total)| Longitud_Total | |||
|---|---|---|---|
| Predictors | Estimates | CI | p |
| (Intercept) | 20.50 | 20.05 – 20.95 | <0.001 |
| Random Effects | |||
| σ2 | 0.21 | ||
| τ00 Id_Random_Foto | 9.15 | ||
| ICC | 0.98 | ||
| N Id_Random_Foto | 175 | ||
| Observations | 525 | ||
| Marginal R2 / Conditional R2 | 0.000 / 0.978 | ||
####Altura_Cola
mixtoalt_cola<-lmer(Altura_Cola~1+(1|Id_Random_Foto),data=datos_morfo)
summary(mixtoalt_cola)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Altura_Cola ~ 1 + (1 | Id_Random_Foto)
## Data: datos_morfo
##
## REML criterion at convergence: -1035.3
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.5269 -0.4111 -0.0115 0.3671 4.7982
##
## Random effects:
## Groups Name Variance Std.Dev.
## Id_Random_Foto (Intercept) 0.057823 0.24046
## Residual 0.001732 0.04161
## Number of obs: 525, groups: Id_Random_Foto, 175
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 1.30511 0.01827 173.99999 71.44 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1icc(mixtoalt_cola, by_group = TRUE)## # ICC by Group
##
## Group | ICC
## ----------------------
## Id_Random_Foto | 0.971tab_model(mixtoalt_cola)| Altura_Cola | |||
|---|---|---|---|
| Predictors | Estimates | CI | p |
| (Intercept) | 1.31 | 1.27 – 1.34 | <0.001 |
| Random Effects | |||
| σ2 | 0.00 | ||
| τ00 Id_Random_Foto | 0.06 | ||
| ICC | 0.97 | ||
| N Id_Random_Foto | 175 | ||
| Observations | 525 | ||
| Marginal R2 / Conditional R2 | 0.000 / 0.971 | ||
####Longitud_Cola
mixtolong_cola<-lmer(Longitud_Cola~1+(1|Id_Random_Foto),data=datos_morfo)
summary(mixtolong_cola)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Longitud_Cola ~ 1 + (1 | Id_Random_Foto)
## Data: datos_morfo
##
## REML criterion at convergence: 1237
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.3111 -0.1864 0.0006 0.1569 5.7480
##
## Random effects:
## Groups Name Variance Std.Dev.
## Id_Random_Foto (Intercept) 4.734 2.1757
## Residual 0.128 0.3577
## Number of obs: 525, groups: Id_Random_Foto, 175
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 14.3668 0.1652 174.0000 86.96 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1icc(mixtolong_cola)## # Intraclass Correlation Coefficient
##
## Adjusted ICC: 0.974
## Unadjusted ICC: 0.974tab_model(mixtolong_cola)| Longitud_Cola | |||
|---|---|---|---|
| Predictors | Estimates | CI | p |
| (Intercept) | 14.37 | 14.04 – 14.69 | <0.001 |
| Random Effects | |||
| σ2 | 0.13 | ||
| τ00 Id_Random_Foto | 4.73 | ||
| ICC | 0.97 | ||
| N Id_Random_Foto | 175 | ||
| Observations | 525 | ||
| Marginal R2 / Conditional R2 | 0.000 / 0.974 | ||
####Perímetro_Musculo
mixtoperim_mus<-lmer(Perim_Musculo1~1+(1|Id_Random_Foto),data=datos_morfo)
summary(mixtoperim_mus)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Perim_Musculo1 ~ 1 + (1 | Id_Random_Foto)
## Data: datos_morfo
##
## REML criterion at convergence: 2094.4
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.9709 -0.2947 -0.0076 0.3260 4.6020
##
## Random effects:
## Groups Name Variance Std.Dev.
## Id_Random_Foto (Intercept) 22.5095 4.7444
## Residual 0.6826 0.8262
## Number of obs: 525, groups: Id_Random_Foto, 175
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 30.6241 0.3605 174.0000 84.96 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1icc(mixtoperim_mus, by_group = TRUE)## # ICC by Group
##
## Group | ICC
## ----------------------
## Id_Random_Foto | 0.971tab_model(mixtoperim_mus)| Perim_Musculo1 | |||
|---|---|---|---|
| Predictors | Estimates | CI | p |
| (Intercept) | 30.62 | 29.92 – 31.33 | <0.001 |
| Random Effects | |||
| σ2 | 0.68 | ||
| τ00 Id_Random_Foto | 22.51 | ||
| ICC | 0.97 | ||
| N Id_Random_Foto | 175 | ||
| Observations | 525 | ||
| Marginal R2 / Conditional R2 | 0.000 / 0.971 | ||
Dada la repetibilidad encontrada en las medidas para todas las variables, usé la mediana de cada individuo y con ese valor realicé el resto de los análisis.
Use la mediana porque esta no es afectada por los valores extremos y se recomienda mas cuando los valores no tienen distribución normal.
###Figuras exploratorias
library(ggplot2)
library(viridis)
library(viridisLite)datos_morfo_med<-read.table("morfo_med_corr.txt",h=T)
head(datos_morfo_med)## Id_Random_Foto TREATMENT Postura Tra_Pos Caja Individuo Longitud_Total Peso
## 1 23 LW C1 LWC1 1 1 24.24000 0.036
## 2 112 LW C1 LWC1 1 2 20.99000 0.027
## 3 27 LW C1 LWC1 1 3 24.96900 0.043
## 4 53 LW C1 LWC1 1 4 23.77800 0.025
## 5 28 LW C1 LWC1 1 5 22.34067 0.036
## 6 162 LW C1 LWC1 1 8 21.71933 0.027
## Altura_Cabeza Altura_Cola Longitud_Cola Perim_Musculo1 Longitud_Cabeza
## 1 2.147000 1.695000 17.33133 37.49933 6.612667
## 2 2.520000 1.315000 14.07900 30.67567 7.106333
## 3 2.454667 1.475333 16.94933 35.16267 8.026667
## 4 2.148000 1.758333 17.08367 37.15100 6.708000
## 5 2.054000 1.493667 15.58667 34.36867 6.583667
## 6 2.071000 1.134667 15.72333 32.51500 6.031667
## L_Longitud_Total L_Peso L_Altura_Cabeza L_Altura_Cola L_Longitud_Cola
## 1 1.384533 -1.443697 0.3318320 0.2291697 1.238832
## 2 1.322012 -1.568636 0.4014005 0.1189258 1.148572
## 3 1.397401 -1.366532 0.3899925 0.1688902 1.229153
## 4 1.376175 -1.602060 0.3320343 0.2451012 1.232581
## 5 1.349096 -1.443697 0.3126004 0.1742537 1.192753
## 6 1.336846 -1.568636 0.3161801 0.0548683 1.196545
## L_Perimetro_Musculo L_Longitud_Cabeza Dim.1 Dim.2 Pcrot
## 1 1.574024 0.8203766 2.7424 -0.6554 1.2360673
## 2 1.486794 0.8516456 1.3287 1.2499 0.5988692
## 3 1.546082 0.9045352 3.2396 0.2805 1.4601928
## 4 1.569971 0.8265931 2.4339 -0.2507 1.0970335
## 5 1.536163 0.8184678 1.7652 -0.4652 0.7956281
## 6 1.512084 0.7804373 0.5242 -0.6478 0.2362660####Altura Cabeza
#para volver una variable categorica:
datos_morfo_med$TREATMENT <- as.factor(datos_morfo_med$TREATMENT)#para cambiar el orden de los niveles tanto en las figuras como en los analisis y que aparezca primero el tratamiento Agua Estable (AE)
datos_morfo_med$TREATMENT <-relevel(datos_morfo_med$TREATMENT, ref="DD")datos_morfo_med$TREATMENT <- as.factor(datos_morfo_med$TREATMENT)
ggplot(datos_morfo_med, aes(x=Postura, y=Altura_Cabeza, fill=TREATMENT)) +
geom_boxplot()+scale_fill_manual(values = c("#B8DE29FF", "#33638DFF","#FDE725FF")) + theme_bw()+labs(x = "CLUTCH",y = "HEAD HEIGHT (mm)")
#los datos:
datos_morfo_trans <- read.table("morfo_med_corr.txt",h=T)library(sjPlot)
library(lme4)
library(lmerTest)#para cambiar el orden de los niveles tanto en las figuras como en los analisis y que aparezca primero el tratamiento Agua Estable (AE)
datos_morfo_trans$TREATMENT<- factor(datos_morfo_trans$TREATMENT , levels = c("HW","DD", "LW"))alt_cabeza<-lmer(L_Altura_Cabeza~TREATMENT+(1|Caja)+(1|Postura),data=datos_morfo_trans)
summary(alt_cabeza)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: L_Altura_Cabeza ~ TREATMENT + (1 | Caja) + (1 | Postura)
## Data: datos_morfo_trans
##
## REML criterion at convergence: -482.8
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.2039 -0.7526 -0.0887 0.6547 2.2227
##
## Random effects:
## Groups Name Variance Std.Dev.
## Caja (Intercept) 0.0001624 0.01274
## Postura (Intercept) 0.0007464 0.02732
## Residual 0.0029887 0.05467
## Number of obs: 175, groups: Caja, 21; Postura, 7
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 0.29423 0.01341 10.39092 21.948 4.81e-10 ***
## TREATMENTDD 0.02663 0.01206 11.84105 2.208 0.0477 *
## TREATMENTLW -0.01429 0.01234 12.89341 -1.158 0.2679
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) TREATMENTD
## TREATMENTDD -0.452
## TREATMENTLW -0.442 0.491tab_model(alt_cabeza)| L_Altura_Cabeza | |||
|---|---|---|---|
| Predictors | Estimates | CI | p |
| (Intercept) | 0.29 | 0.27 – 0.32 | <0.001 |
| TREATMENT [DD] | 0.03 | 0.00 – 0.05 | 0.029 |
| TREATMENT [LW] | -0.01 | -0.04 – 0.01 | 0.249 |
| Random Effects | |||
| σ2 | 0.00 | ||
| τ00 Caja | 0.00 | ||
| τ00 Postura | 0.00 | ||
| ICC | 0.23 | ||
| N Caja | 21 | ||
| N Postura | 7 | ||
| Observations | 175 | ||
| Marginal R2 / Conditional R2 | 0.069 / 0.286 | ||
anova(alt_cabeza)## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## TREATMENT 0.034388 0.017194 2 12.467 5.753 0.01698 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1library(emmeans)
trat.emm.s1 <- emmeans(alt_cabeza, "TREATMENT")
pairs(trat.emm.s1,adjust="fdr")## contrast estimate SE df t.ratio p.value
## HW - DD -0.0266 0.0121 11.3 -2.208 0.0730
## HW - LW 0.0143 0.0124 12.3 1.157 0.2693
## DD - LW 0.0409 0.0123 12.2 3.321 0.0179
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: fdr method for 3 testslibrary(emmeans)
trat.emm.s <- emmeans(alt_cabeza, "TREATMENT")
pairs(trat.emm.s)## contrast estimate SE df t.ratio p.value
## HW - DD -0.0266 0.0121 11.3 -2.208 0.1124
## HW - LW 0.0143 0.0124 12.3 1.157 0.4991
## DD - LW 0.0409 0.0123 12.2 3.321 0.0153
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: tukey method for comparing a family of 3 estimateslibrary(emmeans)
trat.emm.s2 <- emmeans(alt_cabeza, "TREATMENT")
pairs(trat.emm.s2,adjust="none")## contrast estimate SE df t.ratio p.value
## HW - DD -0.0266 0.0121 11.3 -2.208 0.0487
## HW - LW 0.0143 0.0124 12.3 1.157 0.2693
## DD - LW 0.0409 0.0123 12.2 3.321 0.0060
##
## Degrees-of-freedom method: kenward-rogerlibrary(emmeans)
trat.emm.s2 <- emmeans(alt_cabeza, "TREATMENT")
pairs(trat.emm.s2,adjust="bonferroni")## contrast estimate SE df t.ratio p.value
## HW - DD -0.0266 0.0121 11.3 -2.208 0.1461
## HW - LW 0.0143 0.0124 12.3 1.157 0.8079
## DD - LW 0.0409 0.0123 12.2 3.321 0.0179
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: bonferroni method for 3 testslibrary(emmeans)
trat.emm.s2 <- emmeans(alt_cabeza, "TREATMENT")
pairs(trat.emm.s2,adjust="scheffe")## contrast estimate SE df t.ratio p.value
## HW - DD -0.0266 0.0121 11.3 -2.208 0.1317
## HW - LW 0.0143 0.0124 12.3 1.157 0.5299
## DD - LW 0.0409 0.0123 12.2 3.321 0.0197
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: scheffe method with rank 2Comparisons plot:
The blue bars on the plot are the confidence intervals. The red arrowed lines represent a scheme to determine homogeneous groups. If the red lines overlap for two groups, then they are not signficantly different using the method chosen.
The ‘adjust’ argument can take one of several useful methods. ‘tukey’ is default, but others including ‘sidak,’ ‘bonferroni,’ etc can be specified. Specifying ‘none’ produces unadjusted p-values. See help with ‘?emmeans::summary.emmGrid’ for details. Here is an example using the ‘holm’ method of adjustment.
plot(trat.emm.s, comparisons = TRUE)
#para cambiar el orden de los niveles tanto en las figuras como en los analisis y que aparezca primero el tratamiento Agua Estable (AE)
datos_morfo_trans$TREATMENT<- factor(datos_morfo_trans$TREATMENT , levels = c("HW","DD", "LW"))ggplot(datos_morfo_trans, aes(x=TREATMENT, y=L_Altura_Cabeza, fill=TREATMENT)) +
geom_boxplot()+scale_fill_manual(values = c("#33638DFF","#B8DE29FF","#FDE725FF"))+ theme_bw()+labs(x = "Treatment",y = "Head height (mm)")
####Altura Cola
#para cambiar el orden de los niveles tanto en las figuras como en los analisis y que aparezca primero el tratamiento Agua Estable (AE)
datos_morfo_med$TREATMENT <-relevel(datos_morfo_med$TREATMENT, ref="DD")
datos_morfo_trans$TREATMENT<- factor(datos_morfo_trans$TREATMENT , levels = c("HW","DD", "LW"))ggplot(datos_morfo_med, aes(x=Postura, y=Altura_Cola, fill=TREATMENT)) +
geom_boxplot()+scale_fill_manual(values = c("#B8DE29FF", "#33638DFF","#FDE725FF"))+ theme_bw()+labs(x = "Clutch",y = "Tail height (mm)")
ggplot(datos_morfo_trans, aes(x=TREATMENT, y=L_Altura_Cola, fill=TREATMENT)) +
geom_boxplot()+scale_fill_manual(values = c("#33638DFF","#B8DE29FF","#FDE725FF"))+ theme_bw()+labs(x = "Treatment",y = "Head height (mm)")
alt_cola<-lmer(L_Altura_Cola~TREATMENT+(1|Caja)+(1|Postura),data=datos_morfo_trans)
summary(alt_cola)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: L_Altura_Cola ~ TREATMENT + (1 | Caja) + (1 | Postura)
## Data: datos_morfo_trans
##
## REML criterion at convergence: -434.6
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.2606 -0.7026 -0.0209 0.6921 3.8128
##
## Random effects:
## Groups Name Variance Std.Dev.
## Caja (Intercept) 0.0004348 0.02085
## Postura (Intercept) 0.0015837 0.03980
## Residual 0.0038128 0.06175
## Number of obs: 175, groups: Caja, 21; Postura, 7
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 0.09794 0.01876 9.78996 5.220 0.000418 ***
## TREATMENTDD 0.04109 0.01583 11.80604 2.596 0.023652 *
## TREATMENTLW -0.01300 0.01612 12.61909 -0.807 0.434887
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) TREATMENTD
## TREATMENTDD -0.424
## TREATMENTLW -0.416 0.493tab_model(alt_cola)| L_Altura_Cola | |||
|---|---|---|---|
| Predictors | Estimates | CI | p |
| (Intercept) | 0.10 | 0.06 – 0.13 | <0.001 |
| TREATMENT [DD] | 0.04 | 0.01 – 0.07 | 0.010 |
| TREATMENT [LW] | -0.01 | -0.04 – 0.02 | 0.421 |
| Random Effects | |||
| σ2 | 0.00 | ||
| τ00 Caja | 0.00 | ||
| τ00 Postura | 0.00 | ||
| ICC | 0.35 | ||
| N Caja | 21 | ||
| N Postura | 7 | ||
| Observations | 175 | ||
| Marginal R2 / Conditional R2 | 0.084 / 0.401 | ||
anova(alt_cola)## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## TREATMENT 0.047641 0.02382 2 12.296 6.2476 0.01342 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1library(emmeans)
trat.emm <- emmeans(alt_cola, "TREATMENT")
pairs(trat.emm,adjust="fdr")## contrast estimate SE df t.ratio p.value
## HW - DD -0.0411 0.0158 11.5 -2.596 0.0362
## HW - LW 0.0130 0.0161 12.3 0.806 0.4356
## DD - LW 0.0541 0.0161 12.2 3.361 0.0167
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: fdr method for 3 testslibrary(emmeans)
trat.emm <- emmeans(alt_cola, "TREATMENT")
pairs(trat.emm)## contrast estimate SE df t.ratio p.value
## HW - DD -0.0411 0.0158 11.5 -2.596 0.0582
## HW - LW 0.0130 0.0161 12.3 0.806 0.7064
## DD - LW 0.0541 0.0161 12.2 3.361 0.0143
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: tukey method for comparing a family of 3 estimates####Longitud total
ggplot(datos_morfo_med, aes(x=Postura, y=Longitud_Total, fill=TREATMENT)) +
geom_boxplot()+scale_fill_manual(values = c("#B8DE29FF", "#33638DFF","#FDE725FF")) + theme_bw()+labs(x = "Clutch",y = "Total Length (mm)")
ggplot(datos_morfo_trans, aes(x=TREATMENT, y=L_Longitud_Total, fill=TREATMENT)) +
geom_boxplot()+scale_fill_manual(values = c("#33638DFF","#B8DE29FF","#FDE725FF"))+ theme_bw()+labs(x = "Treatment",y = "Total length (mm)")
long_total<-lmer(L_Longitud_Total~TREATMENT+(1|Caja)+(1|Postura),data=datos_morfo_trans)
summary(long_total)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: L_Longitud_Total ~ TREATMENT + (1 | Caja) + (1 | Postura)
## Data: datos_morfo_trans
##
## REML criterion at convergence: -529.6
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.88342 -0.59095 -0.04545 0.60283 2.63211
##
## Random effects:
## Groups Name Variance Std.Dev.
## Caja (Intercept) 0.0006574 0.02564
## Postura (Intercept) 0.0008713 0.02952
## Residual 0.0020930 0.04575
## Number of obs: 175, groups: Caja, 21; Postura, 7
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 1.29852 0.01592 12.01040 81.577 <2e-16 ***
## TREATMENTDD 0.03579 0.01604 11.91700 2.232 0.0456 *
## TREATMENTLW -0.01199 0.01620 12.38584 -0.740 0.4730
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) TREATMENTD
## TREATMENTDD -0.505
## TREATMENTLW -0.500 0.496tab_model(long_total)| L_Longitud_Total | |||
|---|---|---|---|
| Predictors | Estimates | CI | p |
| (Intercept) | 1.30 | 1.27 – 1.33 | <0.001 |
| TREATMENT [DD] | 0.04 | 0.00 – 0.07 | 0.027 |
| TREATMENT [LW] | -0.01 | -0.04 – 0.02 | 0.460 |
| Random Effects | |||
| σ2 | 0.00 | ||
| τ00 Caja | 0.00 | ||
| τ00 Postura | 0.00 | ||
| ICC | 0.42 | ||
| N Caja | 21 | ||
| N Postura | 7 | ||
| Observations | 175 | ||
| Marginal R2 / Conditional R2 | 0.103 / 0.482 | ||
anova(long_total)## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## TREATMENT 0.0199 0.0099501 2 12.206 4.754 0.02973 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1library(emmeans)
trat.emm <- emmeans(long_total, "TREATMENT")
pairs(trat.emm,adjust="fdr")## contrast estimate SE df t.ratio p.value
## HW - DD -0.0358 0.0160 11.7 -2.232 0.0690
## HW - LW 0.0120 0.0162 12.2 0.740 0.4733
## DD - LW 0.0478 0.0162 12.1 2.952 0.0360
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: fdr method for 3 testslibrary(emmeans)
trat.emm <- emmeans(long_total, "TREATMENT")
pairs(trat.emm)## contrast estimate SE df t.ratio p.value
## HW - DD -0.0358 0.0160 11.7 -2.232 0.1068
## HW - LW 0.0120 0.0162 12.2 0.740 0.7451
## DD - LW 0.0478 0.0162 12.1 2.952 0.0300
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: tukey method for comparing a family of 3 estimates####Longitud Cola
datos_morfo_med$TREATMENT <-relevel(datos_morfo_med$TREATMENT, ref="DD")
ggplot(datos_morfo_med, aes(x=Postura, y=Longitud_Cola, fill=TREATMENT)) +
geom_boxplot()+scale_fill_manual(values = c("#B8DE29FF", "#33638DFF","#FDE725FF")) + theme_bw()+labs(x = "Clutch",y = "Tail Length (mm)")
ggplot(datos_morfo_trans, aes(x=TREATMENT, y=L_Longitud_Cola, fill=TREATMENT)) +
geom_boxplot()+scale_fill_manual(values = c("#33638DFF","#B8DE29FF","#FDE725FF"))+ theme_bw()+labs(x = "Treatment",y = "Tail length (mm)")
long_cola<-lmer(L_Longitud_Cola~TREATMENT+(1|Caja)+(1|Postura),data=datos_morfo_trans)
summary(long_cola)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: L_Longitud_Cola ~ TREATMENT + (1 | Caja) + (1 | Postura)
## Data: datos_morfo_trans
##
## REML criterion at convergence: -518.7
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.14532 -0.61867 0.02129 0.58558 2.50560
##
## Random effects:
## Groups Name Variance Std.Dev.
## Caja (Intercept) 0.0007150 0.02674
## Postura (Intercept) 0.0008673 0.02945
## Residual 0.0022315 0.04724
## Number of obs: 175, groups: Caja, 21; Postura, 7
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 1.14365 0.01623 12.34127 70.476 <2e-16 ***
## TREATMENTDD 0.03760 0.01668 11.92947 2.254 0.0438 *
## TREATMENTLW -0.01307 0.01685 12.39238 -0.776 0.4523
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) TREATMENTD
## TREATMENTDD -0.515
## TREATMENTLW -0.510 0.496tab_model(long_cola)| L_Longitud_Cola | |||
|---|---|---|---|
| Predictors | Estimates | CI | p |
| (Intercept) | 1.14 | 1.11 – 1.18 | <0.001 |
| TREATMENT [DD] | 0.04 | 0.00 – 0.07 | 0.025 |
| TREATMENT [LW] | -0.01 | -0.05 – 0.02 | 0.439 |
| Random Effects | |||
| σ2 | 0.00 | ||
| τ00 Caja | 0.00 | ||
| τ00 Postura | 0.00 | ||
| ICC | 0.41 | ||
| N Caja | 21 | ||
| N Postura | 7 | ||
| Observations | 175 | ||
| Marginal R2 / Conditional R2 | 0.109 / 0.479 | ||
anova(long_cola)## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## TREATMENT 0.021964 0.010982 2 12.214 4.9213 0.02707 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1library(emmeans)
trat.emm <- emmeans(long_cola, "TREATMENT")
pairs(trat.emm,adjust="fdr")## contrast estimate SE df t.ratio p.value
## HW - DD -0.0376 0.0167 11.7 -2.254 0.0663
## HW - LW 0.0131 0.0169 12.2 0.776 0.4527
## DD - LW 0.0507 0.0168 12.1 3.010 0.0323
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: fdr method for 3 testslibrary(emmeans)
trat.emm <- emmeans(long_cola, "TREATMENT")
pairs(trat.emm)## contrast estimate SE df t.ratio p.value
## HW - DD -0.0376 0.0167 11.7 -2.254 0.1029
## HW - LW 0.0131 0.0169 12.2 0.776 0.7242
## DD - LW 0.0507 0.0168 12.1 3.010 0.0270
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: tukey method for comparing a family of 3 estimates####Perímetro Músculo
ggplot(datos_morfo_med, aes(x=Postura, y=L_Perimetro_Musculo, fill=TREATMENT)) +
geom_boxplot()+scale_fill_manual(values = c("#B8DE29FF", "#33638DFF","#FDE725FF")) + theme_bw()+labs(x = "Clutch",y = "Muscle perimeter (mm)")
ggplot(datos_morfo_trans, aes(x=TREATMENT, y=L_Perimetro_Musculo, fill=TREATMENT)) +
geom_boxplot()+scale_fill_manual(values = c("#33638DFF","#B8DE29FF","#FDE725FF"))+ theme_bw()+labs(x = "Treatment",y = "Muscle perimeter (mm)")
peri_mus<-lmer(L_Perimetro_Musculo~TREATMENT+(1|Caja)+(1|Postura),data=datos_morfo_trans)
summary(peri_mus)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: L_Perimetro_Musculo ~ TREATMENT + (1 | Caja) + (1 | Postura)
## Data: datos_morfo_trans
##
## REML criterion at convergence: -508.8
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.2480 -0.5780 0.0161 0.6813 2.5950
##
## Random effects:
## Groups Name Variance Std.Dev.
## Caja (Intercept) 0.0007589 0.02755
## Postura (Intercept) 0.0008978 0.02996
## Residual 0.0023654 0.04864
## Number of obs: 175, groups: Caja, 21; Postura, 7
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 1.47229 0.01662 12.43165 88.589 <2e-16 ***
## TREATMENTDD 0.03721 0.01718 11.91715 2.166 0.0513 .
## TREATMENTLW -0.01337 0.01735 12.37915 -0.770 0.4555
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) TREATMENTD
## TREATMENTDD -0.518
## TREATMENTLW -0.513 0.496tab_model(peri_mus)| L_Perimetro_Musculo | |||
|---|---|---|---|
| Predictors | Estimates | CI | p |
| (Intercept) | 1.47 | 1.44 – 1.51 | <0.001 |
| TREATMENT [DD] | 0.04 | 0.00 – 0.07 | 0.032 |
| TREATMENT [LW] | -0.01 | -0.05 – 0.02 | 0.442 |
| Random Effects | |||
| σ2 | 0.00 | ||
| τ00 Caja | 0.00 | ||
| τ00 Postura | 0.00 | ||
| ICC | 0.41 | ||
| N Caja | 21 | ||
| N Postura | 7 | ||
| Observations | 175 | ||
| Marginal R2 / Conditional R2 | 0.103 / 0.473 | ||
anova(peri_mus)## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## TREATMENT 0.021776 0.010888 2 12.201 4.6029 0.0324 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1library(emmeans)
trat.emm <- emmeans(peri_mus, "TREATMENT")
pairs(trat.emm,adjust="fdr")## contrast estimate SE df t.ratio p.value
## HW - DD -0.0372 0.0172 11.7 -2.166 0.0776
## HW - LW 0.0134 0.0174 12.2 0.770 0.4558
## DD - LW 0.0506 0.0173 12.1 2.917 0.0384
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: fdr method for 3 testslibrary(emmeans)
trat.emm <- emmeans(peri_mus, "TREATMENT")
pairs(trat.emm)## contrast estimate SE df t.ratio p.value
## HW - DD -0.0372 0.0172 11.7 -2.166 0.1191
## HW - LW 0.0134 0.0174 12.2 0.770 0.7274
## DD - LW 0.0506 0.0173 12.1 2.917 0.0319
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: tukey method for comparing a family of 3 estimates####Longitud_Cabeza
ggplot(datos_morfo_med, aes(x=Postura, y=Longitud_Cabeza, fill=TREATMENT)) +
geom_boxplot()+scale_fill_manual(values = c("#B8DE29FF", "#33638DFF","#FDE725FF")) + theme_bw()+labs(x = "Clutch",y = "Head length (mm)")
ggplot(datos_morfo_trans, aes(x=TREATMENT, y=L_Longitud_Cabeza, fill=TREATMENT)) +
geom_boxplot()+scale_fill_manual(values = c("#33638DFF","#B8DE29FF","#FDE725FF"))+ theme_bw()+labs(x = "Treatment",y = "Head Length (mm)")
long_cabeza<-lmer(L_Longitud_Cabeza~TREATMENT+(1|Caja)+(1|Postura),data=datos_morfo_trans)
summary(long_cabeza)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: L_Longitud_Cabeza ~ TREATMENT + (1 | Caja) + (1 | Postura)
## Data: datos_morfo_trans
##
## REML criterion at convergence: -489.9
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.19812 -0.64890 -0.05108 0.54024 2.79254
##
## Random effects:
## Groups Name Variance Std.Dev.
## Caja (Intercept) 0.0005131 0.02265
## Postura (Intercept) 0.0010190 0.03192
## Residual 0.0027145 0.05210
## Number of obs: 175, groups: Caja, 21; Postura, 7
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 0.77541 0.01625 11.08876 47.703 3.46e-14 ***
## TREATMENTDD 0.03239 0.01538 11.87412 2.106 0.0572 .
## TREATMENTLW -0.01431 0.01560 12.51147 -0.917 0.3763
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) TREATMENTD
## TREATMENTDD -0.475
## TREATMENTLW -0.468 0.495tab_model(long_cabeza)| L_Longitud_Cabeza | |||
|---|---|---|---|
| Predictors | Estimates | CI | p |
| (Intercept) | 0.78 | 0.74 – 0.81 | <0.001 |
| TREATMENT [DD] | 0.03 | 0.00 – 0.06 | 0.037 |
| TREATMENT [LW] | -0.01 | -0.05 – 0.02 | 0.360 |
| Random Effects | |||
| σ2 | 0.00 | ||
| τ00 Caja | 0.00 | ||
| τ00 Postura | 0.00 | ||
| ICC | 0.36 | ||
| N Caja | 21 | ||
| N Postura | 7 | ||
| Observations | 175 | ||
| Marginal R2 / Conditional R2 | 0.083 / 0.414 | ||
anova(long_cabeza)## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## TREATMENT 0.02583 0.012915 2 12.262 4.7578 0.02955 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1library(emmeans)
trat.emm <- emmeans(long_cabeza, "TREATMENT")
pairs(trat.emm,adjust="fdr")## contrast estimate SE df t.ratio p.value
## HW - DD -0.0324 0.0154 11.6 -2.106 0.0866
## HW - LW 0.0143 0.0156 12.2 0.917 0.3769
## DD - LW 0.0467 0.0156 12.1 2.997 0.0330
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: fdr method for 3 testslibrary(emmeans)
trat.emm <- emmeans(long_cabeza, "TREATMENT")
pairs(trat.emm)## contrast estimate SE df t.ratio p.value
## HW - DD -0.0324 0.0154 11.6 -2.106 0.1317
## HW - LW 0.0143 0.0156 12.2 0.917 0.6402
## DD - LW 0.0467 0.0156 12.1 2.997 0.0276
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: tukey method for comparing a family of 3 estimates####Peso
ggplot(datos_morfo_med, aes(x=Postura, y=Peso, fill=TREATMENT)) +
geom_boxplot()+scale_fill_manual(values = c("#B8DE29FF", "#33638DFF","#FDE725FF"))+ theme_bw()+labs(x = "CLUTCH",y = "BODY WEIGHT (g)")
ggplot(datos_morfo_trans, aes(x=TREATMENT, y=L_Peso, fill=TREATMENT)) +
geom_boxplot()+scale_fill_manual(values = c("#33638DFF","#B8DE29FF","#FDE725FF"))+ theme_bw()+labs(x = "Treatment",y = "Body Weight (g)")
peso<-lmer(L_Peso~TREATMENT+(1|Caja)+(1|Postura),data=datos_morfo_trans)
summary(peso)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: L_Peso ~ TREATMENT + (1 | Caja) + (1 | Postura)
## Data: datos_morfo_trans
##
## REML criterion at convergence: -211.1
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.7125 -0.4795 0.0077 0.5846 3.4941
##
## Random effects:
## Groups Name Variance Std.Dev.
## Caja (Intercept) 0.008624 0.09287
## Postura (Intercept) 0.009873 0.09936
## Residual 0.012490 0.11176
## Number of obs: 175, groups: Caja, 21; Postura, 7
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -1.68116 0.05340 12.01296 -31.484 6.49e-13 ***
## TREATMENTDD 0.10379 0.05365 11.94032 1.935 0.0771 .
## TREATMENTLW -0.04935 0.05395 12.20180 -0.915 0.3780
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) TREATMENTD
## TREATMENTDD -0.503
## TREATMENTLW -0.500 0.498tab_model(peso)| L_Peso | |||
|---|---|---|---|
| Predictors | Estimates | CI | p |
| (Intercept) | -1.68 | -1.79 – -1.58 | <0.001 |
| TREATMENT [DD] | 0.10 | -0.00 – 0.21 | 0.055 |
| TREATMENT [LW] | -0.05 | -0.16 – 0.06 | 0.362 |
| Random Effects | |||
| σ2 | 0.01 | ||
| τ00 Caja | 0.01 | ||
| τ00 Postura | 0.01 | ||
| ICC | 0.60 | ||
| N Caja | 21 | ||
| N Postura | 7 | ||
| Observations | 175 | ||
| Marginal R2 / Conditional R2 | 0.117 / 0.644 | ||
anova(peso)## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## TREATMENT 0.10531 0.052654 2 12.103 4.2158 0.0408 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1library(emmeans)
trat.emm <- emmeans(peso, "TREATMENT")
pairs(trat.emm,adjust="fdr")## contrast estimate SE df t.ratio p.value
## HW - DD -0.1038 0.0536 11.8 -1.935 0.1159
## HW - LW 0.0494 0.0540 12.1 0.915 0.3782
## DD - LW 0.1531 0.0539 12.1 2.840 0.0445
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: fdr method for 3 testslibrary(emmeans)
trat.emm <- emmeans(peso, "TREATMENT")
pairs(trat.emm)## contrast estimate SE df t.ratio p.value
## HW - DD -0.1038 0.0536 11.8 -1.935 0.1720
## HW - LW 0.0494 0.0540 12.1 0.915 0.6416
## DD - LW 0.1531 0.0539 12.1 2.840 0.0367
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: tukey method for comparing a family of 3 estimates###Test de normalidad y Homogenidad de varianzas
library(tidyverse)
library(rstatix)
library(ggpubr)Normalidad:
#####Altura Cabeza
#Shapiro_test
datos_morfo_med %>%group_by(TREATMENT) %>%
shapiro_test(Altura_Cabeza)## # A tibble: 3 × 4
## TREATMENT variable statistic p
## <fct> <chr> <dbl> <dbl>
## 1 DD Altura_Cabeza 0.959 0.0411
## 2 HW Altura_Cabeza 0.964 0.0762
## 3 LW Altura_Cabeza 0.989 0.898#####Longitud_Total
#Shapiro_test
datos_morfo_med %>%group_by(TREATMENT) %>%
shapiro_test(Longitud_Total)## # A tibble: 3 × 4
## TREATMENT variable statistic p
## <fct> <chr> <dbl> <dbl>
## 1 DD Longitud_Total 0.866 0.00000828
## 2 HW Longitud_Total 0.958 0.0371
## 3 LW Longitud_Total 0.987 0.834#####Altura_Cola
#Shapiro_test
datos_morfo_med %>%group_by(TREATMENT) %>%
shapiro_test(Altura_Cola)## # A tibble: 3 × 4
## TREATMENT variable statistic p
## <fct> <chr> <dbl> <dbl>
## 1 DD Altura_Cola 0.946 0.00904
## 2 HW Altura_Cola 0.978 0.352
## 3 LW Altura_Cola 0.962 0.0890#####Longitud_Cola
#Shapiro_test
datos_morfo_med %>%group_by(TREATMENT) %>%
shapiro_test(Longitud_Cola)## # A tibble: 3 × 4
## TREATMENT variable statistic p
## <fct> <chr> <dbl> <dbl>
## 1 DD Longitud_Cola 0.885 0.0000344
## 2 HW Longitud_Cola 0.981 0.477
## 3 LW Longitud_Cola 0.988 0.843#####Perímetro_Musculo
#Shapiro_test
datos_morfo_med %>%group_by(TREATMENT) %>%
shapiro_test(Perim_Musculo1)## # A tibble: 3 × 4
## TREATMENT variable statistic p
## <fct> <chr> <dbl> <dbl>
## 1 DD Perim_Musculo1 0.888 0.0000429
## 2 HW Perim_Musculo1 0.975 0.248
## 3 LW Perim_Musculo1 0.992 0.975#####Longitud_Cabeza
#Shapiro_test
datos_morfo_med %>%group_by(TREATMENT) %>%
shapiro_test(Longitud_Cabeza)## # A tibble: 3 × 4
## TREATMENT variable statistic p
## <fct> <chr> <dbl> <dbl>
## 1 DD Longitud_Cabeza 0.922 0.000843
## 2 HW Longitud_Cabeza 0.899 0.000123
## 3 LW Longitud_Cabeza 0.943 0.0121#####Peso
#Shapiro_test
datos_morfo_med %>%group_by(TREATMENT) %>%
shapiro_test(Peso)## # A tibble: 3 × 4
## TREATMENT variable statistic p
## <fct> <chr> <dbl> <dbl>
## 1 DD Peso 0.757 0.0000000110
## 2 HW Peso 0.939 0.00496
## 3 LW Peso 0.907 0.000497Normalidad variables transformadas:
#los datos:
datos_morfo_med<- read.table("morfo_med_corr.txt",h=T)#####log Altura Cabeza
#Shapiro_test
datos_morfo_med %>%group_by(TREATMENT) %>%
shapiro_test(L_Altura_Cabeza)## # A tibble: 3 × 4
## TREATMENT variable statistic p
## <chr> <chr> <dbl> <dbl>
## 1 DD L_Altura_Cabeza 0.985 0.676
## 2 HW L_Altura_Cabeza 0.977 0.308
## 3 LW L_Altura_Cabeza 0.986 0.760#####Log Longitud_Total
#Shapiro_test
datos_morfo_med %>%group_by(TREATMENT) %>%
shapiro_test(L_Longitud_Total)## # A tibble: 3 × 4
## TREATMENT variable statistic p
## <chr> <chr> <dbl> <dbl>
## 1 DD L_Longitud_Total 0.911 0.000306
## 2 HW L_Longitud_Total 0.977 0.319
## 3 LW L_Longitud_Total 0.989 0.891#####Log Altura_Cola
#Shapiro_test
datos_morfo_med %>%group_by(TREATMENT) %>%
shapiro_test(L_Altura_Cola)## # A tibble: 3 × 4
## TREATMENT variable statistic p
## <chr> <chr> <dbl> <dbl>
## 1 DD L_Altura_Cola 0.972 0.180
## 2 HW L_Altura_Cola 0.992 0.965
## 3 LW L_Altura_Cola 0.976 0.343#####Log Longitud_Cola
#Shapiro_test
datos_morfo_med %>%group_by(TREATMENT) %>%
shapiro_test(L_Longitud_Cola)## # A tibble: 3 × 4
## TREATMENT variable statistic p
## <chr> <chr> <dbl> <dbl>
## 1 DD L_Longitud_Cola 0.926 0.00122
## 2 HW L_Longitud_Cola 0.990 0.887
## 3 LW L_Longitud_Cola 0.984 0.661#####log Perímetro_Musculo
#Shapiro_test
datos_morfo_med %>%group_by(TREATMENT) %>%
shapiro_test(L_Perimetro_Musculo)## # A tibble: 3 × 4
## TREATMENT variable statistic p
## <chr> <chr> <dbl> <dbl>
## 1 DD L_Perimetro_Musculo 0.934 0.00265
## 2 HW L_Perimetro_Musculo 0.987 0.795
## 3 LW L_Perimetro_Musculo 0.984 0.665#####Log Longitud_Cabeza
#Shapiro_test
datos_morfo_med%>%group_by(TREATMENT) %>%
shapiro_test(L_Longitud_Cabeza)## # A tibble: 3 × 4
## TREATMENT variable statistic p
## <chr> <chr> <dbl> <dbl>
## 1 DD L_Longitud_Cabeza 0.954 0.0212
## 2 HW L_Longitud_Cabeza 0.946 0.00990
## 3 LW L_Longitud_Cabeza 0.966 0.130#####Log Peso
#Shapiro_test
datos_morfo_med %>%group_by(TREATMENT) %>%
shapiro_test(L_Peso)## # A tibble: 3 × 4
## TREATMENT variable statistic p
## <chr> <chr> <dbl> <dbl>
## 1 DD L_Peso 0.939 0.00425
## 2 HW L_Peso 0.972 0.188
## 3 LW L_Peso 0.977 0.383####Homogenidad de varianza (de cada variable por separado):
library(car)Levene’s test: A robust alternative to the Bartlett’s test that is less sensitive to departures from normality.
Fligner-Killeen’s test: a non-parametric test which is very robust against departures from normality.
#1. Las variables independientes deben ser categoricas
datos_morfo_med$TREATMENT <- as.factor(datos_morfo_med$TREATMENT)
datos_morfo_med$Postura<- as.factor(datos_morfo_med$Postura)
datos_morfo_med$Caja<- as.factor(datos_morfo_med$Caja)Longitud_Total
# Levene's test with multiple independent variables
leveneTest( Longitud_Total~TREATMENT*Postura*Caja, data = datos_morfo_med)## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 20 0.5524 0.9385
## 154Longitud_Cola
# Levene's test with multiple independent variables
leveneTest( Longitud_Cola~TREATMENT, data = datos_morfo_med)## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 3.5681 0.03032 *
## 172
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Longitud_Cabeza
# Levene's test with multiple independent variables
leveneTest( Longitud_Cabeza~TREATMENT, data = datos_morfo_med)## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 0.6404 0.5283
## 172Altura_Cabeza
# Levene's test with multiple independent variables
leveneTest( Altura_Cabeza~TREATMENT, data = datos_morfo_med)## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 1.583 0.2083
## 172Altura_Cola
# Levene's test with multiple independent variables
leveneTest( Altura_Cola~TREATMENT, data = datos_morfo_med)## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 3.529 0.03148 *
## 172
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Perímetro_Músculo
# Levene's test with multiple independent variables
leveneTest(Perim_Musculo1~TREATMENT, data = datos_morfo_med)## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 3.2137 0.04263 *
## 172
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Peso
# Levene's test with multiple independent variables
leveneTest(Peso~TREATMENT, data = datos_morfo_med)## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 5.4672 0.004989 **
## 172
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1####Homogenidad de varianza variables transformadas (de cada variable por separado):
library(car)Levene’s test: A robust alternative to the Bartlett’s test that is less sensitive to departures from normality.
Fligner-Killeen’s test: a non-parametric test which is very robust against departures from normality.
#1. Las variables independientes deben ser categoricas
datos_morfo_med$TREATMENT <- as.factor(datos_morfo_med$TREATMENT)
datos_morfo_med$Postura<- as.factor(datos_morfo_med$Postura)
datos_morfo_med$Caja<- as.factor(datos_morfo_med$Caja)log(Longitud_Total)
# Levene's test with multiple independent variables
leveneTest( L_Longitud_Total~TREATMENT*Postura*Caja, data = datos_morfo_med)## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 20 0.6051 0.9047
## 154log Longitud_Cola
# Levene's test with multiple independent variables
leveneTest( L_Longitud_Cola~TREATMENT, data = datos_morfo_med)## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 1.8937 0.1536
## 172log Longitud_Cabeza
# Levene's test with multiple independent variables
leveneTest( L_Longitud_Cabeza~TREATMENT, data = datos_morfo_med)## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 0.176 0.8388
## 172log Altura_Cabeza
# Levene's test with multiple independent variables
leveneTest( L_Altura_Cabeza~TREATMENT, data = datos_morfo_med)## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 0.6533 0.5216
## 172log Altura_Cola
# Levene's test with multiple independent variables
leveneTest( L_Altura_Cola~TREATMENT, data = datos_morfo_med)## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 1.6916 0.1873
## 172log Perímetro_Músculo
# Levene's test with multiple independent variables
leveneTest(L_Perimetro_Musculo~TREATMENT, data = datos_morfo_med)## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 1.7292 0.1805
## 172log Peso
# Levene's test with multiple independent variables
leveneTest(L_Peso~TREATMENT, data = datos_morfo_med)## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 2 1.7396 0.1787
## 172###Figuras exploratorias datos transformados
library(ggplot2)#los datos:
datos_morfo_trans <- read.table("morfo_med_corr.txt",h=T)Altura_cabeza
#para volver una variable categorica:
datos_morfo_trans$TREATMENT <- as.factor(datos_morfo_trans$TREATMENT)#para cambiar el orden de los niveles tanto en las figuras como en los analisis y que aparezca primero el tratamiento DRY DOWN
datos_morfo_trans$TREATMENT <-relevel(datos_morfo_trans$TREATMENT, ref="DD")ggplot(datos_morfo_trans, aes(x=Postura, y=L_Altura_Cabeza, fill=TREATMENT)) + geom_boxplot()+scale_fill_manual(values = c("#B8DE29FF", "#33638DFF","#FDE725FF")) + theme_bw()+labs(x = "Clutch",y = "log(Head height) (mm)")
Altura_cola
ggplot(datos_morfo_trans, aes(x=Postura, y=L_Altura_Cola, fill=TREATMENT)) + geom_boxplot()+scale_fill_manual(values = c("#B8DE29FF", "#33638DFF","#FDE725FF")) + theme_bw()+labs(x = "Clutch",y = "log(Tail height) (mm)")
Longitud total
ggplot(datos_morfo_trans, aes(x=Postura, y=L_Longitud_Total, fill=TREATMENT)) + geom_boxplot()+scale_fill_manual(values = c("#B8DE29FF", "#33638DFF","#FDE725FF")) + theme_bw()+labs(x = "Clutch",y = "log(Total length) (mm)")
Longitud cola
ggplot(datos_morfo_trans, aes(x=Postura, y=L_Longitud_Cola, fill=TREATMENT)) + geom_boxplot()+scale_fill_manual(values = c("#B8DE29FF", "#33638DFF","#FDE725FF")) + theme_bw()+labs(x = "Clutch",y = "log(Tail length) (mm)")
Perímetro músculo
ggplot(datos_morfo_trans, aes(x=Postura, y=L_Perimetro_Musculo, fill=TREATMENT)) + geom_boxplot()+scale_fill_manual(values = c("#B8DE29FF", "#33638DFF","#FDE725FF")) + theme_bw()+labs(x = "Clutch",y = "log(Muscle perimeter) (mm)")
Longitud Cabeza
ggplot(datos_morfo_trans, aes(x=Postura, y=L_Longitud_Cabeza, fill=TREATMENT)) + geom_boxplot()+scale_fill_manual(values = c("#B8DE29FF", "#33638DFF","#FDE725FF")) + theme_bw()+labs(x = "Clutch",y = "log(Head length) (mm)")
Peso
ggplot(datos_morfo_trans, aes(x=Postura, y=L_Peso, fill=TREATMENT)) + geom_boxplot()+scale_fill_manual(values = c("#B8DE29FF", "#33638DFF","#FDE725FF")) + theme_bw()+labs(x = "Clutch",y = "log(Weight) (g)")
Análisis de componentes pricipales variables transformadas
#los datos:
datos_morfo_trans <- read.table("morfo_med_corr.txt",h=T)####Correlación entre variables transformadas
morfo<-datos_morfo_trans[,15:20] #primero los individuos activos y luego las variables
head(morfo)## L_Peso L_Altura_Cabeza L_Altura_Cola L_Longitud_Cola L_Perimetro_Musculo
## 1 -1.443697 0.3318320 0.2291697 1.238832 1.574024
## 2 -1.568636 0.4014005 0.1189258 1.148572 1.486794
## 3 -1.366532 0.3899925 0.1688902 1.229153 1.546082
## 4 -1.602060 0.3320343 0.2451012 1.232581 1.569971
## 5 -1.443697 0.3126004 0.1742537 1.192753 1.536163
## 6 -1.568636 0.3161801 0.0548683 1.196545 1.512084
## L_Longitud_Cabeza
## 1 0.8203766
## 2 0.8516456
## 3 0.9045352
## 4 0.8265931
## 5 0.8184678
## 6 0.7804373#####La Correlación
corrl<-round(cor(morfo),3) Las figuras de las correlaciones con coeficientes
library(corrplot)
library(ggcorrplot)
ggcorrplot(corrl, hc.order = TRUE, type = "lower",
lab = TRUE, colors = c("#E3E418FF","#440154FF","#75D054FF"))
####El PCA con las variables transformadas
#####Los paquetes:
para hacer el test de barlett y el analisis de componentes: library(psych)
para graficar que variables pesan mas en que PC library(corrplot)
for ggplot2-based visualization: library(factoextra)
library(factoextra)
library(FactoMineR)
library(psych)#####El test de bartlett: Ho:the sample correlation matrix came from a multivariate normal population in which the variables of interest are independent.
Rejection of the hypothesis is taken as an indication that the data are appropriate for analysis.
(CHARLES D. DZIUBAN, 1974) Bartlett’s test of sphericity tests the hypothesis that your correlation matrix is an identity matrix,which would indicate that your variables are unrelated and therefore unsuitable for structure detection. Small values (less than 0.05) of the significance level indicate that a factor analysis may be useful with your data.
cortest.bartlett(corrl,n=175)## $chisq
## [1] 1364.241
##
## $p.value
## [1] 8.14692e-282
##
## $df
## [1] 15KMO: Test alternativo al test de barlett. Measure of Sampling Adequacy
This is just a function of the squared elements of the ‘image’ matrix compared to the squares of the original correlations. The overall MSA as well as estimates for each item are found.The index is known as the Kaiser-Meyer-Olkin (KMO) index.Kaiser-Meyer-Olkin (KMO). KMO Test is a measure of how suited your data is for Factor Analysis.
KMO(morfo)## Kaiser-Meyer-Olkin factor adequacy
## Call: KMO(r = morfo)
## Overall MSA = 0.88
## MSA for each item =
## L_Peso L_Altura_Cabeza L_Altura_Cola L_Longitud_Cola
## 0.96 0.92 0.95 0.79
## L_Perimetro_Musculo L_Longitud_Cabeza
## 0.79 0.93#####Análisis de Componentes principales:
#PCA usando factorMineR y factorextra
elpca<-PCA(morfo, scale.unit = TRUE,ncp=6,graph = FALSE)
summary(elpca)##
## Call:
## PCA(X = morfo, scale.unit = TRUE, ncp = 6, graph = FALSE)
##
##
## Eigenvalues
## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5 Dim.6
## Variance 4.894 0.373 0.294 0.240 0.184 0.015
## % of var. 81.569 6.217 4.901 3.997 3.063 0.252
## Cumulative % of var. 81.569 87.786 92.687 96.685 99.748 100.000
##
## Individuals (the 10 first)
## Dist Dim.1 ctr cos2 Dim.2 ctr cos2
## 1 | 2.914 | 2.742 0.878 0.886 | -0.655 0.658 0.051 |
## 2 | 2.001 | 1.329 0.206 0.441 | 1.250 2.393 0.390 |
## 3 | 3.369 | 3.240 1.225 0.925 | 0.281 0.121 0.007 |
## 4 | 2.735 | 2.434 0.692 0.792 | -0.251 0.096 0.008 |
## 5 | 1.919 | 1.765 0.364 0.846 | -0.465 0.331 0.059 |
## 6 | 1.240 | 0.524 0.032 0.179 | -0.648 0.643 0.273 |
## 7 | 1.159 | 0.642 0.048 0.307 | 0.187 0.054 0.026 |
## 8 | 1.992 | 1.964 0.450 0.973 | 0.205 0.064 0.011 |
## 9 | 2.701 | 2.448 0.700 0.821 | 0.393 0.236 0.021 |
## 10 | 0.844 | -0.685 0.055 0.658 | 0.331 0.168 0.154 |
## Dim.3 ctr cos2
## 1 0.648 0.816 0.049 |
## 2 -0.538 0.562 0.072 |
## 3 -0.416 0.337 0.015 |
## 4 0.542 0.571 0.039 |
## 5 0.455 0.403 0.056 |
## 6 -0.700 0.952 0.319 |
## 7 -0.433 0.365 0.140 |
## 8 0.061 0.007 0.001 |
## 9 -0.033 0.002 0.000 |
## 10 0.182 0.065 0.047 |
##
## Variables
## Dim.1 ctr cos2 Dim.2 ctr cos2 Dim.3
## L_Peso | 0.890 16.172 0.791 | -0.174 8.141 0.030 | 0.162
## L_Altura_Cabeza | 0.871 15.509 0.759 | 0.391 40.917 0.153 | -0.132
## L_Altura_Cola | 0.872 15.553 0.761 | 0.141 5.361 0.020 | 0.427
## L_Longitud_Cola | 0.940 18.053 0.884 | -0.269 19.461 0.073 | -0.114
## L_Perimetro_Musculo | 0.947 18.319 0.897 | -0.243 15.883 0.059 | -0.118
## L_Longitud_Cabeza | 0.896 16.394 0.802 | 0.195 10.237 0.038 | -0.204
## ctr cos2
## L_Peso 8.876 0.026 |
## L_Altura_Cabeza 5.947 0.017 |
## L_Altura_Cola 61.946 0.182 |
## L_Longitud_Cola 4.408 0.013 |
## L_Perimetro_Musculo 4.700 0.014 |
## L_Longitud_Cabeza 14.124 0.042 |#diagrama de codos
fviz_eig(elpca, addlabels = TRUE, ylim = c(0, 90))
#proporcion de varianza explicada:
eig.val <- get_eigenvalue(elpca)variables <- get_pca_var(elpca)# Contribucion de las variables en los componentes
corrplot(variables$contrib, is.corr=FALSE)
#el biplot de las cargas
fviz_pca_var (elpca, col.var = "contrib",
gradient.cols=c("#481567FF","#287D8EFF","#3CBB75FF","#B8DE29FF"))
#Para exportar los componentes:
elpca$ind$coord## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
## 1 2.74236350 -0.6554173613 0.6480155952 -0.329143448 -0.078599094
## 2 1.32866320 1.2498788200 -0.5375808491 0.459301061 0.416956899
## 3 3.23961296 0.2805463680 -0.4161732621 0.758537271 0.017174702
## 4 2.43390052 -0.2506912864 0.5420189102 -0.950385692 -0.537703528
## 5 1.76519634 -0.4651881116 0.4554623958 0.340299584 -0.025592661
## 6 0.52418441 -0.6478130762 -0.7000372933 0.152446666 0.557450618
## 7 0.64228706 0.1870495739 -0.4334845788 0.764352301 -0.306586799
## 8 1.96406254 0.2047403915 0.0608398683 -0.109563001 -0.185070328
## 9 2.44828107 0.3926745510 -0.0334100949 1.065121941 -0.100174494
## 10 -0.68455411 0.3313160092 0.1822646967 0.097068877 -0.284832083
## 11 -1.79771241 0.8630487372 -0.0759901156 -0.954599389 0.519807556
## 12 -1.73113028 -0.1559480094 0.7576616044 -0.815249641 -0.704788696
## 13 0.22695506 0.6367779312 -0.6253380680 -0.667374957 -0.905298396
## 14 -4.83371905 1.5488360930 0.2713597367 0.366059536 0.292284162
## 15 -1.62018582 0.1366680436 0.4951317579 -0.054150769 0.252956499
## 16 1.31881684 0.8378695228 -0.6600172811 -0.286361726 -0.626040565
## 17 -1.64794651 0.2312191650 0.0928712000 -0.042313860 -0.117659375
## 18 -2.67496873 0.0895723063 0.7208582309 -0.305916428 0.012046576
## 19 1.55982400 0.3414877222 -0.7468740874 -0.336412413 -0.446610043
## 20 -3.88270835 1.0168925745 -0.3734532783 -0.187213059 0.441850872
## 21 -2.88304939 0.1333168671 0.3835939052 -0.498190017 -0.099566510
## 22 1.35366703 0.4926308547 -0.8200329949 -0.578453606 -0.966176930
## 23 -2.33021431 -0.1019342410 -0.9629324352 -0.217400237 -0.901238307
## 24 -2.83050067 0.4748869985 0.2525912724 -0.602389482 -0.162882009
## 25 -1.63991518 -0.2371396508 -1.0480573526 -0.243295529 0.378810207
## 26 -0.43327776 0.4110483230 -0.6106493642 0.007690087 -0.420683063
## 27 -2.24373600 -0.6798671878 1.5750364323 -0.721384807 -0.454144598
## 28 -1.79035695 0.7524407009 0.7677303285 0.119054034 -0.017527568
## 29 -4.00379965 -0.0129118854 0.5352778061 0.180583617 -0.431337995
## 30 -0.90506368 0.6371444910 -0.2876781702 0.225792426 0.297837122
## 31 -1.53057224 -0.2622948928 0.0926983268 0.066400163 0.289455046
## 32 -0.30200141 0.0873254389 0.7503350071 -0.142322706 0.642232880
## 33 -3.60152118 -0.3940982660 0.0670718496 -0.521572376 -0.010686727
## 34 -0.45847044 -0.0290749505 0.4747039358 0.138948004 0.126835745
## 35 -0.92800498 0.1040522202 1.2005184720 -0.130690250 -0.233043216
## 36 -0.72105547 0.6566377207 -0.1027221396 0.035936551 0.669238702
## 37 -2.11713235 -0.0139208859 -0.0533228102 0.240646793 0.533248754
## 38 -4.33564030 0.3901250551 0.2320327691 0.187747956 0.441904636
## 39 1.00203704 -0.5359453571 0.1102655284 0.111678229 0.037554733
## 40 -0.25715891 -0.1132322226 0.3299238503 -0.362091359 0.331303902
## 41 0.85166847 -0.4420184135 0.1936991497 -0.294496490 0.746243859
## 42 -0.26026300 -0.1202505715 -0.4039922906 0.123393971 0.568848944
## 43 0.83715520 -0.2164505970 0.3074413018 -0.386940273 0.285766198
## 44 -1.56352318 0.5393511916 -0.2157140757 0.556668939 -0.106581637
## 45 -2.24995704 -0.6475364454 0.4691218729 -0.062133611 -0.391710813
## 46 -0.62995719 -0.4575022433 -0.2811759968 -0.074786632 0.088216769
## 47 0.01844112 0.0213951725 0.1920855710 -0.570212765 -0.349167714
## 48 0.06428149 -0.7797749136 -0.7308613387 0.243824935 0.789776832
## 49 -2.94113751 -1.0902693622 -0.5988439965 0.115830801 -0.359053551
## 50 -0.21935540 -0.4870750081 -0.4926308388 -0.003097605 -0.025728205
## 51 -0.97481372 -0.7866013235 -0.6899091266 -0.533916933 0.343820619
## 52 -1.35565003 -0.9457291633 -0.5036848250 0.769887664 -0.322753705
## 53 -1.60371556 -0.6555117424 -0.4214481959 0.399965256 -0.151732356
## 54 -3.54087886 -0.7462846386 0.2131653031 0.785163082 -0.507780988
## 55 7.02725419 -0.4973493240 0.2847739294 0.339050686 0.069577791
## 56 6.51775575 -0.3322648927 -0.0548539021 0.375194224 0.120427291
## 57 6.04633322 0.5509026557 -0.2029868865 -0.554611579 0.283229509
## 58 5.92958589 -1.4849034363 0.4049040602 0.548923161 -0.233370930
## 59 5.87492976 -0.5259243123 0.5934506464 0.430612284 -0.399057261
## 60 5.81156300 -1.1262717809 -0.6027587904 0.425920381 0.164815919
## 61 1.91414409 -0.8276371499 0.2205054510 -0.013943982 0.038361339
## 62 5.16089712 -0.2051753369 0.2531205562 0.492090187 0.366304979
## 63 1.33671666 0.6409999128 0.5547494526 0.069730984 0.070858022
## 64 0.88272900 0.3535639471 -0.4551862895 0.812740425 -0.518119681
## 65 1.66520239 0.1496060602 -0.0628009658 0.069684991 -0.157465678
## 66 4.10001058 0.2468816616 -0.0920036122 -0.981184167 -0.868564029
## 67 1.42046341 -0.0875966768 -0.1769806778 -0.034412997 -0.028460226
## 68 4.21590231 0.0217200507 0.3817159240 -0.291533538 -0.308516470
## 69 0.19596442 0.1350209470 0.5120999193 0.662451171 0.255039303
## 70 0.96248651 -0.3353938284 1.2246146825 0.245200835 -0.458500305
## 71 -1.59836012 0.1207419905 0.5922396617 0.347548679 -0.269743138
## 72 0.70473782 0.5499556000 -0.1168560858 -0.144630134 0.333811320
## 73 0.74796979 0.0754943263 -0.3426685698 0.014859733 -0.192194998
## 74 1.36080966 0.5110409776 -0.5348099486 -0.574167080 -0.816871028
## 75 -1.32983479 0.3205504523 -0.0143589579 0.271304782 0.208965882
## 76 -1.90863595 0.2378658986 -0.4424559593 -0.599027121 -0.215952242
## 77 -0.28360886 0.7775195741 -0.1537350892 0.362598721 0.029007223
## 78 -0.39941777 0.8434665670 -0.5516728787 -0.274338367 -0.441582926
## 79 0.05269538 0.1856479430 0.2765780358 0.276364401 -0.076379937
## 80 -0.68225004 0.5413195622 0.2786041835 1.039432305 -0.481219515
## 81 -1.77378948 0.1327766946 -1.3471314121 -0.607326010 -0.271479090
## 82 1.05676240 1.2207026676 -0.4924615385 -0.882791651 -0.177924425
## 83 -0.27525291 0.0995153608 -1.7220896756 -0.486254695 0.180897467
## 84 -1.22032824 -0.3291631103 0.1154420957 -0.062292923 -0.207713059
## 85 -2.34174274 -0.5731779026 0.1401932460 0.235422909 0.208827915
## 86 1.31350300 0.3021591147 -0.0158878125 -0.883887448 -0.743818751
## 87 -1.98365587 -0.4109234699 -0.6294907803 0.413376010 0.219678397
## 88 -3.02530497 -0.2370094652 0.1160098247 0.251677610 -0.001890055
## 89 -2.53217744 0.0169389602 -0.3268087312 0.521463688 -0.495258475
## 90 1.09366998 0.7182226482 1.3818615048 0.394498926 -0.015354072
## 91 1.35021369 0.5697812242 0.4667729758 0.580224380 0.546076476
## 92 0.04505815 -0.3141835614 0.0579144506 0.102680884 0.315486797
## 93 -0.22550603 -0.2601319229 0.3909579667 0.315193382 0.231301096
## 94 4.15038941 0.5298250017 -0.4111213802 -0.534537248 0.362853800
## 95 -1.23731660 -0.2740744016 0.1277258363 0.111016258 0.191598502
## 96 1.59923209 0.6285490377 0.1381031651 0.752349818 0.074512957
## 97 -1.28895539 0.9490028300 0.1125500214 1.250249476 -0.081469513
## 98 0.26769125 1.1582607543 -0.3317988388 0.919261685 -0.055083764
## 99 3.04564775 0.7793979195 0.5075706578 -0.568477820 0.872275076
## 100 0.42260148 -0.1817777928 0.3647068256 -0.856322471 0.472433263
## 101 0.67843488 0.3562208159 0.2431131659 -0.383433126 0.921653247
## 102 2.89968971 -0.1933147967 -0.0051745360 -0.235757267 0.437700555
## 103 -0.02121621 -0.1518774325 -0.0923869769 0.143565239 0.280582875
## 104 0.20952036 -1.0201756175 0.2800104662 -0.302226149 -0.324949750
## 105 0.01028680 -0.6629385541 -0.0427910682 0.306345757 0.041203692
## 106 -1.84291074 -0.6429867669 -0.6239474394 -0.175313739 0.249829960
## 107 0.90630068 -0.2816024777 0.8376566123 -0.590418815 0.138585715
## 108 0.70397768 -0.6024263406 -0.2817654442 -0.209250325 -0.001344548
## 109 2.08298533 -0.2416942660 0.5420374411 0.106615724 -0.001148914
## 110 1.23698627 -0.3480166944 0.2199694243 0.360195446 -0.197087211
## 111 1.59959761 -1.9003667027 -0.6303880993 -0.829535940 1.773014515
## 112 -2.07658367 -0.6247838489 -0.3161417719 -0.182601084 0.279224304
## 113 1.07403374 -0.6606277412 0.4594269401 -0.158005743 -0.113439401
## 114 -0.25441690 -1.1535926917 0.3358863017 0.053165270 -0.650212262
## 115 -2.61787392 -0.6677463277 0.3839579640 -0.618861166 -0.748822764
## 116 0.16644592 -0.2704737029 0.3130186708 0.308266233 -0.367655952
## 117 -0.17882789 -0.4277775877 -0.1851728697 -0.435202816 -0.088771097
## 118 1.42839858 0.1426621510 0.4928804849 0.097361084 0.043726111
## 119 3.37179575 0.9166981827 -1.0730218832 -0.738827799 -1.220447826
## 120 0.10645466 -0.8328074211 -0.4969140603 0.457301304 -0.376525435
## 121 -0.61554080 -0.3166518640 -0.1290885005 0.388425741 0.077459531
## 122 0.01581282 0.4335553187 0.7476480324 0.005732329 -0.268582550
## 123 1.20636402 0.6185887406 0.4937120631 -0.336998593 0.535893376
## 124 -3.30889866 0.3854305042 0.7151904462 0.023710852 -0.195323035
## 125 1.47634075 -0.2742630006 -0.3911487357 -0.239593228 -0.164086875
## 126 1.79696776 -0.2617583234 -0.0003603035 -0.173639715 0.250923067
## 127 0.59863724 0.0671864752 -0.1346416473 0.816197394 -0.010438765
## 128 3.26317644 0.2258902683 -0.4886353711 -0.381703410 -0.154855342
## 129 1.37899789 -0.8199253892 0.0502805148 0.049989179 0.021791434
## 130 -1.49705175 0.3299138318 0.0147933576 0.085755233 0.006687474
## 131 1.13849382 0.1576861543 0.6751348050 0.551815215 -0.222793721
## 132 1.89201654 0.7679483540 0.4757330475 -0.151177390 0.719032633
## 133 -0.31568792 0.1506749843 1.0731293570 -0.167172893 -0.058576317
## 134 2.48539825 -0.4660107388 0.0823674596 0.597172975 0.018690471
## 135 -0.63565019 -0.4671184611 -0.3185478591 0.244576359 0.405649127
## 136 -0.62685721 0.8119930245 -0.0010919830 0.436827655 0.452810857
## 137 -1.52528509 -0.3388334473 0.1105707656 0.225970008 0.242492460
## 138 -2.74515560 -0.2587877422 0.1049615026 -0.332008026 -0.126125106
## 139 2.21549312 1.4922232508 -1.0520946325 -0.193781201 -0.622862681
## 140 -1.64005340 -0.1245504602 -0.1451445950 0.456558316 0.288713324
## 141 -1.14628746 0.1651978083 0.4823181671 -0.125366945 -0.078193638
## 142 1.04944565 0.5193304723 0.3343118215 0.236547115 0.369009630
## 143 -3.34367592 0.3018754601 -0.5278815336 0.773229990 0.249609169
## 144 -1.45122127 0.5632998622 0.5547870880 -0.682262027 -0.082774557
## 145 0.99806353 -0.0721809677 -0.1972799018 -0.331647452 0.520023000
## 146 -0.96468916 -0.1199569687 0.0945364147 -0.987872003 -0.169255968
## 147 -1.48605633 -0.7943983772 -0.5945279365 -0.660233467 -0.203931677
## 148 -3.07283675 0.8842737062 1.0572091462 -0.162768819 -0.498844364
## 149 -0.46127422 1.6037972203 -0.9465286421 -0.724060646 0.429457989
## 150 -2.93948821 0.0083507071 0.5812699078 -0.088296300 -0.072649795
## 151 2.24890774 0.0906313842 -0.2863385679 -0.521211822 0.239445232
## 152 4.54290561 0.4057355602 -0.5234267263 -0.069725999 -0.411841308
## 153 0.20465024 -0.1576837846 -0.2845711201 0.457893099 0.111398642
## 154 -0.77743673 -0.7298231097 -0.1024249002 0.026487715 0.316839809
## 155 2.02516058 0.4685589146 0.1202141394 1.256470097 -0.396731879
## 156 -3.36484589 -0.4159558999 0.0667821197 0.672251809 -0.124146942
## 157 -0.41152061 0.5657279126 0.1600053608 0.532205302 0.307754554
## 158 0.08767643 -0.0010629522 0.4300202820 0.555012375 -0.142469940
## 159 1.88849957 0.5929659614 -0.3617385341 1.083547516 0.220211339
## 160 -0.60117896 -0.9987519203 0.3837405933 0.032054322 -0.453675001
## 161 -1.28286525 0.4618197618 -0.7593069819 -0.374226881 0.547531217
## 162 -1.55182835 0.2099055661 -0.2002879036 -0.230590854 0.739043761
## 163 -0.64333371 1.3838835803 -0.5745688503 1.063393503 -0.453411228
## 164 -0.17796956 0.4091591500 0.5115079820 -0.547622773 0.843362890
## 165 -1.21178942 -0.3527341311 0.3893850684 0.282174875 -0.105642762
## 166 -2.25245452 -0.8743855021 0.3242485657 -0.409690145 -0.646033045
## 167 -1.94306834 0.0007871925 0.2813213524 0.115262937 -0.030359984
## 168 -2.24352243 -0.0469419787 -0.5173338383 0.022535107 0.240385818
## 169 -3.80975027 -0.1499940585 -0.9139855945 -0.131669796 0.981185875
## 170 0.91763019 -0.0306949201 0.1834965648 -0.364964118 0.889360054
## 171 1.69438597 0.7940780944 1.3942833145 -1.328676335 -0.048898411
## 172 -2.04722471 -0.5980788563 -0.9808553289 -0.089850998 -0.243818410
## 173 -1.29885142 -0.8238108798 -0.1061243617 -0.175052492 -0.104675729
## 174 -0.54575631 -1.3630059411 -0.7964102122 -0.119115832 -0.190206300
## 175 -1.63352957 -1.2284047372 -0.8828761081 0.524076624 -0.769894854
## Dim.6
## 1 -0.069804195
## 2 -0.065629977
## 3 0.168417502
## 4 -0.083687885
## 5 -0.164653675
## 6 0.136048674
## 7 0.171540579
## 8 0.130169510
## 9 -0.059282068
## 10 -0.102276691
## 11 -0.099429054
## 12 0.113446316
## 13 -0.209080176
## 14 0.089758727
## 15 0.104715849
## 16 0.236738593
## 17 -0.139721852
## 18 0.052193884
## 19 0.162273770
## 20 0.047634703
## 21 0.087143212
## 22 -0.151390867
## 23 -0.138052875
## 24 -0.061615759
## 25 -0.126836541
## 26 -0.269304076
## 27 0.008054668
## 28 -0.007127540
## 29 0.043307163
## 30 -0.119542713
## 31 -0.110871180
## 32 -0.015810604
## 33 -0.161692447
## 34 -0.170933896
## 35 -0.147100059
## 36 0.189205830
## 37 0.087724929
## 38 0.083922027
## 39 -0.121225581
## 40 0.011460515
## 41 0.047719168
## 42 0.045506067
## 43 -0.067078175
## 44 0.086008973
## 45 0.060919586
## 46 0.089128669
## 47 -0.077476844
## 48 -0.066201432
## 49 -0.050031996
## 50 -0.037952538
## 51 -0.087479014
## 52 0.108509275
## 53 0.118431109
## 54 0.115996188
## 55 -0.104857055
## 56 -0.143051132
## 57 -0.057843951
## 58 -0.089944520
## 59 -0.060708371
## 60 0.148925338
## 61 0.112578972
## 62 0.127312866
## 63 -0.147610086
## 64 -0.229436712
## 65 -0.025553910
## 66 0.093157309
## 67 -0.190140463
## 68 0.002774389
## 69 0.113880521
## 70 0.110023059
## 71 0.077297749
## 72 -0.091966087
## 73 -0.148503263
## 74 0.152454310
## 75 -0.140559259
## 76 -0.089176088
## 77 -0.113745730
## 78 0.154186490
## 79 0.118054110
## 80 0.112108553
## 81 0.048374584
## 82 -0.039090944
## 83 -0.034563689
## 84 -0.075380733
## 85 -0.028933849
## 86 0.075534182
## 87 0.107655049
## 88 0.071048700
## 89 0.087929417
## 90 -0.097532385
## 91 0.117892598
## 92 -0.060405739
## 93 -0.059606962
## 94 0.195104316
## 95 -0.037137086
## 96 0.196638582
## 97 0.096362686
## 98 0.119106290
## 99 0.155842918
## 100 0.041829621
## 101 -0.152623392
## 102 0.143751714
## 103 -0.225045596
## 104 0.073466254
## 105 0.120352542
## 106 0.101555252
## 107 0.068497475
## 108 -0.088346430
## 109 -0.222538459
## 110 -0.173163374
## 111 0.150242633
## 112 -0.083337963
## 113 0.147646517
## 114 0.132758219
## 115 0.070101783
## 116 -0.109457183
## 117 0.020866008
## 118 -0.046719363
## 119 0.233162094
## 120 -0.197066784
## 121 -0.129032037
## 122 0.046109002
## 123 0.074148449
## 124 0.018282311
## 125 -0.088871947
## 126 -0.028899579
## 127 -0.116390336
## 128 -0.084626882
## 129 -0.069977675
## 130 -0.081260784
## 131 0.117199695
## 132 0.093238575
## 133 0.018794661
## 134 0.171919374
## 135 -0.129437752
## 136 -0.229511147
## 137 -0.029516199
## 138 -0.201717809
## 139 0.140972833
## 140 0.065537638
## 141 0.171710980
## 142 0.085642920
## 143 -0.015647935
## 144 -0.026728865
## 145 0.269025679
## 146 -0.068584171
## 147 0.043071758
## 148 0.048898174
## 149 0.172188404
## 150 0.055049633
## 151 -0.123925702
## 152 -0.130879214
## 153 -0.093371834
## 154 -0.133110662
## 155 -0.237360399
## 156 -0.123612849
## 157 0.347638960
## 158 0.101553278
## 159 0.125973792
## 160 0.092938398
## 161 -0.112052868
## 162 -0.120154848
## 163 -0.302393324
## 164 -0.151632553
## 165 0.066639954
## 166 0.074670664
## 167 0.055183072
## 168 -0.102824417
## 169 0.043289969
## 170 -0.046121115
## 171 -0.062519846
## 172 0.040901144
## 173 0.111437693
## 174 0.118517574
## 175 0.120515849# Export into a TXT file
write.infile(elpca$ind$coord, "pca_trans.txt", sep = "\t")
# Export into a CSV file
write.infile(elpca$ind$coord, "pca_trans.csv", sep = ";")#####Rotando los componentes
Abdi & Williams, 2010, Principal component analysis: After the number of components has been determined, and in order to facilitate the interpretation,the analysis often involves a rotation of the components that were retained. principal: Standardized loadings (pattern matrix) based upon correlation matrix
morforot <- principal(morfo,1,rotate="varimax")morforot1<-principal(morfo, rotate="varimax", nfactors=1, scores=TRUE)cargas_rot<-morforot1$loadingsPara extraer los componentes ya rotados
componentes_rot<-morforot1$scores[1:175,]
write.csv(componentes_rot,file ="morfo_rot_trans.csv")####Figuras con el PCrotado para ver el efecto del tratamiento
datos_morfo_trans$TREATMENT <- as.factor(datos_morfo_trans$TREATMENT)
datos_morfo_trans$TREATMENT <-relevel(datos_morfo_trans$TREATMENT, ref="DD")
ggplot(datos_morfo_trans, aes(x=Postura, y=Pcrot, fill=TREATMENT)) +
geom_boxplot()+scale_fill_manual(values = c("#B8DE29FF", "#33638DFF","#FDE725FF"))+ theme_bw()+labs(x = "CLUTCH",y = "Rotated PC (Body Size)")
###El modelo mixto para evaluar el efecto del tratamiento de desecación:
Aqui otra forma de presentar los resultados, calculando el R2 del modelo y haciendo un ajuste en la significancia:
Como crear una tabla con los resultados de los análisis de modelos mixtos
*Daniel Lüdecke
The marginal R-squared considers only the variance of the fixed effects, while the conditional R-squared takes both the fixed and random effects into account.
The p-value is a simple approximation, based on the t-statistics and using the normal distribution function. A more precise p-value can be computed using p.val = “kr”. In this case, which only applies to linear mixed models, the computation of p-values is based on conditional F-tests with Kenward-Roger approximation for the degrees of freedom (using the using the pbkrtest-package). Note that here the computation is more time consuming and thus not used as default. You can also display the approximated degrees of freedom with show.df.
Para calcular el R*2 del modelo
Nakagawa S, Johnson P, Schielzeth H (2017) The coefficient of determination R2 and intra-class correlation coefficient from generalized linear mixed-effects models revisted and expanded. J. R. Soc. Interface 14. doi: 10.1098/rsif.2017.0213.
library(sjPlot)
library(lme4)
library(lmerTest)El modelo mixto con PC rotado
datos_morfo_trans$TREATMENT<-as.factor(datos_morfo_trans$TREATMENT)
datos_morfo_trans$TREATMENT <-relevel(datos_morfo_trans$TREATMENT, ref="HW")
mixto1pc1<-lmer(Pcrot~TREATMENT+(1|Caja)+(1|Postura),data=datos_morfo_trans)
summary(mixto1pc1)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Pcrot ~ TREATMENT + (1 | Caja) + (1 | Postura)
## Data: datos_morfo_trans
##
## REML criterion at convergence: 421.9
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.57051 -0.62221 -0.03631 0.61635 2.54916
##
## Random effects:
## Groups Name Variance Std.Dev.
## Caja (Intercept) 0.1713 0.4139
## Postura (Intercept) 0.2663 0.5161
## Residual 0.5249 0.7245
## Number of obs: 175, groups: Caja, 21; Postura, 7
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -0.1370 0.2670 11.3731 -0.513 0.6176
## TREATMENTDD 0.5940 0.2576 11.9251 2.306 0.0399 *
## TREATMENTLW -0.2403 0.2601 12.3826 -0.924 0.3732
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) TREATMENTD
## TREATMENTDD -0.483
## TREATMENTLW -0.479 0.496anova(mixto1pc1)## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## TREATMENT 5.7703 2.8851 2 12.207 5.4971 0.01985 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1tab_model(mixto1pc1)| Pcrot | |||
|---|---|---|---|
| Predictors | Estimates | CI | p |
| (Intercept) | -0.14 | -0.66 – 0.39 | 0.608 |
| TREATMENT [DD] | 0.59 | 0.09 – 1.10 | 0.022 |
| TREATMENT [LW] | -0.24 | -0.75 – 0.27 | 0.357 |
| Random Effects | |||
| σ2 | 0.52 | ||
| τ00 Caja | 0.17 | ||
| τ00 Postura | 0.27 | ||
| ICC | 0.45 | ||
| N Caja | 21 | ||
| N Postura | 7 | ||
| Observations | 175 | ||
| Marginal R2 / Conditional R2 | 0.114 / 0.517 | ||
library(sjPlot)
library(lme4)
library(lmerTest)
tab_model(mixto1pc1)| Pcrot | |||
|---|---|---|---|
| Predictors | Estimates | CI | p |
| (Intercept) | -0.14 | -0.66 – 0.39 | 0.608 |
| TREATMENT [DD] | 0.59 | 0.09 – 1.10 | 0.022 |
| TREATMENT [LW] | -0.24 | -0.75 – 0.27 | 0.357 |
| Random Effects | |||
| σ2 | 0.52 | ||
| τ00 Caja | 0.17 | ||
| τ00 Postura | 0.27 | ||
| ICC | 0.45 | ||
| N Caja | 21 | ||
| N Postura | 7 | ||
| Observations | 175 | ||
| Marginal R2 / Conditional R2 | 0.114 / 0.517 | ||
Haciendo el ajuste en la significancia:
“The p-value is a simple approximation, based on the t-statistics and using the normal distribution function. A more precise p-value can be computed using p.val =”kr”. In this case, which only applies to linear mixed models, the computation of p-values is based on conditional F-tests with Kenward-Roger approximation for the degrees of freedom (using the using the pbkrtest-package). Note that here the computation is more time consuming and thus not used as default. You can also display the approximated degrees of freedom with show.df.”
tab_model(mixto1pc1, p.val = "kr", show.df = TRUE)| Pcrot | ||||
|---|---|---|---|---|
| Predictors | Estimates | CI | p | df |
| (Intercept) | -0.14 | -0.72 – 0.45 | 0.618 | 11.40 |
| TREATMENT [DD] | 0.59 | 0.03 – 1.16 | 0.040 | 11.71 |
| TREATMENT [LW] | -0.24 | -0.81 – 0.33 | 0.374 | 12.16 |
| Random Effects | ||||
| σ2 | 0.52 | |||
| τ00 Caja | 0.17 | |||
| τ00 Postura | 0.27 | |||
| ICC | 0.45 | |||
| N Caja | 21 | |||
| N Postura | 7 | |||
| Observations | 175 | |||
| Marginal R2 / Conditional R2 | 0.114 / 0.517 | |||
ggplot(datos_morfo_trans, aes(x=TREATMENT, y=Pcrot, fill=TREATMENT)) +
geom_boxplot()+scale_fill_manual(values = c("#33638DFF","#B8DE29FF","#FDE725FF"))+ theme_bw()+labs(x = "Treatment",y = "Body Size (Rotated component)")
library(yarrr)datos_morfo_trans$TREATMENT<-as.factor(datos_morfo_trans$TREATMENT)
datos_morfo_trans$TREATMENT <-relevel(datos_morfo_trans$TREATMENT, ref="HW")pirateplot(formula = Pcrot~TREATMENT,
data =datos_morfo_trans,
theme = 1,pal= c("#33638DFF","#B8DE29FF","#FDE725FF"),xlab = "Treatment", ylab="Body Size (PC rotated)",point.o = 1, point.pch = 20, point.col = c("#33638DFF","#B8DE29FF","#FDE725FF"),bean.f.o = 0.1,inf.f.col = c("#33638DFF","#B8DE29FF","#FDE725FF"),inf.f.o =0.9)
library(emmeans)
trat.emm <- emmeans(mixto1pc1, "TREATMENT")
pairs(trat.emm,adjust="fdr")## contrast estimate SE df t.ratio p.value
## HW - DD -0.594 0.258 11.7 -2.306 0.0604
## HW - LW 0.240 0.260 12.2 0.924 0.3736
## DD - LW 0.834 0.260 12.1 3.210 0.0222
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: fdr method for 3 testslibrary(lsmeans)
lsmeans::lsmeans(mixto1pc1, pairwise~TREATMENT, adjust="fdr")## $lsmeans
## TREATMENT lsmean SE df lower.CL upper.CL
## HW -0.137 0.267 11.4 -0.722 0.448
## DD 0.457 0.267 11.3 -0.128 1.042
## LW -0.377 0.269 11.8 -0.965 0.211
##
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df t.ratio p.value
## HW - DD -0.594 0.258 11.7 -2.306 0.0604
## HW - LW 0.240 0.260 12.2 0.924 0.3736
## DD - LW 0.834 0.260 12.1 3.210 0.0222
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: fdr method for 3 testslibrary(emmeans)
emmeans(mixto1pc1, list(pairwise ~ TREATMENT), adjust = "tukey")## $`emmeans of TREATMENT`
## TREATMENT emmean SE df lower.CL upper.CL
## HW -0.137 0.267 11.4 -0.722 0.448
## DD 0.457 0.267 11.3 -0.128 1.042
## LW -0.377 0.269 11.8 -0.965 0.211
##
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $`pairwise differences of TREATMENT`
## 1 estimate SE df t.ratio p.value
## HW - DD -0.594 0.258 11.7 -2.306 0.0943
## HW - LW 0.240 0.260 12.2 0.924 0.6362
## DD - LW 0.834 0.260 12.1 3.210 0.0188
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: tukey method for comparing a family of 3 estimatestrat.emm <- emmeans(mixto1pc1, "TREATMENT")
pairs(trat.emm,adjust="none") #Least significant difference## contrast estimate SE df t.ratio p.value
## HW - DD -0.594 0.258 11.7 -2.306 0.0402
## HW - LW 0.240 0.260 12.2 0.924 0.3736
## DD - LW 0.834 0.260 12.1 3.210 0.0074
##
## Degrees-of-freedom method: kenward-rogerlibrary(emmeans)
trat.emm <- emmeans(mixto1pc1, "TREATMENT")
pairs(trat.emm,adjust="none") #Least significant difference## contrast estimate SE df t.ratio p.value
## HW - DD -0.594 0.258 11.7 -2.306 0.0402
## HW - LW 0.240 0.260 12.2 0.924 0.3736
## DD - LW 0.834 0.260 12.1 3.210 0.0074
##
## Degrees-of-freedom method: kenward-rogerlibrary(emmeans)
trat.emm <- emmeans(mixto1pc1, "TREATMENT")
pairs(trat.emm)## contrast estimate SE df t.ratio p.value
## HW - DD -0.594 0.258 11.7 -2.306 0.0943
## HW - LW 0.240 0.260 12.2 0.924 0.6362
## DD - LW 0.834 0.260 12.1 3.210 0.0188
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: tukey method for comparing a family of 3 estimates1.Prueba de homogeneidad del modelo mixto rotado pc1
mixto1pc1<-lmer(Pcrot~TREATMENT+(1|Caja)+(1|Postura),data=datos_morfo_trans)
par(mfrow=c(2,2)) # init 4 charts in 1 panel
plot(mixto1pc1)
Prueba de homogeneidad sin considerar los factores aleatorios (lineal model)
trat_pcrot<-lm(Pcrot~TREATMENT,data=datos_morfo_trans)
par(mfrow=c(2,2)) # init 4 charts in 1 panel
plot(trat_pcrot)
Prueba de homogeneidad Breusch Pagan:
library(olsrr) #para el analisis de homog. var de breusch_pagan. No funciona si el modelo no es lm??
datos_morfo_trans <- read.table("morfo_med_completa_expcasa_trans.txt",h=T)
#hice una columna combinando tratamiento y postura. En esa variable Tra_Pos quedaron 21 categorias.
tra_pos_pcrot<-lm(Pcrot~Tra_Pos, data=datos_morfo_trans)
plot(fitted(tra_pos_pcrot),residuals(tra_pos_pcrot))
ols_test_breusch_pagan(tra_pos_pcrot)##
## Breusch Pagan Test for Heteroskedasticity
## -----------------------------------------
## Ho: the variance is constant
## Ha: the variance is not constant
##
## Data
## ---------------------------------
## Response : Pcrot
## Variables: fitted values of Pcrot
##
## Test Summary
## -----------------------------
## DF = 1
## Chi2 = 0.01377641
## Prob > Chi2 = 0.9065645par(mfrow=c(2,2)) # init 4 charts in 1 panel
plot(tra_pos_pcrot)
Otra forma de evaluar la homogeneidad propuesta por Michael Palmeri
library(lmerTest)
library(lme4)datos_morfo_trans <- read.table("morfo_med_corr.txt",h=T)
mixto1pc1<-lmer(Pcrot~TREATMENT+(1|Caja)+(1|Postura),data=datos_morfo_trans)#extracts the residuals and places them in a new column in our original data table
datos_morfo_trans$mixto1pc1.Res<-residuals(mixto1pc1)#creates a new column with the absolute value of the residuals
datos_morfo_trans$Abs.mixto1pc1.Res<-abs(datos_morfo_trans$mixto1pc1.Res)#squares the absolute values of the residuals to provide the more robust estimate
datos_morfo_trans$mixto1pc1.Res2 <- datos_morfo_trans$Abs.mixto1pc1.Res^2 #ANOVA of the squared residuals
Levene.mixto1pc1 <- lm(mixto1pc1.Res2 ~Caja+Postura, data=datos_morfo_trans)
anova(Levene.mixto1pc1) #displays the results## Analysis of Variance Table
##
## Response: mixto1pc1.Res2
## Df Sum Sq Mean Sq F value Pr(>F)
## Caja 1 0.049 0.04949 0.1074 0.7435
## Postura 6 2.171 0.36190 0.7855 0.5824
## Residuals 167 76.944 0.460742.Evaluando los supuestos de los modelos mixtos:
evaluando los supuestos en modelos mixtos
Francisco Juretig. “R Statistics Cookbook.”: Mixed models can be impacted greatly by outliers. Even a minor contamination causes major estimation problems. Solución: evaluar que tan robusto es el modelo.
#para volver una variable categorica:
datos_morfo_trans$TREATMENT<-as.factor(datos_morfo_trans$TREATMENT)Re-estimando el modelo con “robust errors”
Based on central model approach: That is, we assume the model to be true, but a part of the data to possibly be contaminated. Irrespective of this contamination, we would like to estimate the parameters that define the central model and these estimates should be only minimally influenced by the contamination.
A robust estimation method should give reasonable results in the presence of such contaminated error distributions, while still maintain- ing efficiency in case there is no contamination.
library(robustlmm)datos_morfo_trans$TREATMENT <-relevel(datos_morfo_trans$TREATMENT, ref="DD")mixto1pc1_robust <-rlmer(Pcrot~TREATMENT+(1|Caja)+(1|Postura),data=datos_morfo_trans)
summary(mixto1pc1_robust)## Robust linear mixed model fit by DAStau
## Formula: Pcrot ~ TREATMENT + (1 | Caja) + (1 | Postura)
## Data: datos_morfo_trans
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.75764 -0.63946 -0.01392 0.56091 2.71445
##
## Random effects:
## Groups Name Variance Std.Dev.
## Caja (Intercept) 0.1009 0.3176
## Postura (Intercept) 0.2182 0.4671
## Residual 0.4980 0.7057
## Number of obs: 175, groups: Caja, 21; Postura, 7
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.3333 0.2378 1.401
## TREATMENTHW -0.5191 0.2184 -2.377
## TREATMENTLW -0.7522 0.2210 -3.404
##
## Correlation of Fixed Effects:
## (Intr) TREATMENTH
## TREATMENTHW -0.458
## TREATMENTLW -0.452 0.493
##
## Robustness weights for the residuals:
## 139 weights are ~= 1. The remaining 36 ones are summarized as
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.488 0.687 0.794 0.788 0.923 0.999
##
## Robustness weights for the random effects:
## 26 weights are ~= 1. The remaining 2 ones are
## 8 22
## 0.407 0.691
##
## Rho functions used for fitting:
## Residuals:
## eff: smoothed Huber (k = 1.345, s = 10)
## sig: smoothed Huber, Proposal 2 (k = 1.345, s = 10)
## Random Effects, variance component 1 (Caja):
## eff: smoothed Huber (k = 1.345, s = 10)
## vcp: smoothed Huber, Proposal 2 (k = 1.345, s = 10)
## Random Effects, variance component 2 (Postura):
## eff: smoothed Huber (k = 1.345, s = 10)
## vcp: smoothed Huber, Proposal 2 (k = 1.345, s = 10)library(sjPlot)
library(lme4)
tab_model(mixto1pc1_robust)| Pcrot | |||
|---|---|---|---|
| Predictors | Estimates | CI | p |
| (Intercept) | 0.33 | -0.13 – 0.80 | 0.161 |
| TREATMENT [HW] | -0.52 | -0.95 – -0.09 | 0.017 |
| TREATMENT [LW] | -0.75 | -1.19 – -0.32 | 0.001 |
| Random Effects | |||
| σ2 | 0.50 | ||
| τ00 Caja | 0.10 | ||
| τ00 Postura | 0.22 | ||
| ICC | 0.39 | ||
| N Caja | 21 | ||
| N Postura | 7 | ||
| Observations | 175 | ||
| Marginal R2 / Conditional R2 | 0.108 / 0.457 | ||
plot(mixto1pc1_robust)
Evaluando los residuos del modelo
Excerpt From: Francisco Juretig. “R Statistics Cookbook.” Apple Books.
“In the fitted-residual relationship; ideally, there should be no structure there. We then plotted boxplots for the residuals in each group.
Excerpt From: Francisco Juretig. “R Statistics Cookbook.” Apple Books.
Efecto aleatorio:La postura
plot(mixto1pc1, resid(., scaled=TRUE) ~ fitted(.) |Postura, abline = 0)
#The following output shows the scaled residuals for each group:
plot(mixto1pc1, Postura ~ resid(., scaled=TRUE)) 
“Now, let’s plot the fitted versus the actuals. There should be a positive linear relationship here. If this weren’t the case, it would imply that we were missing some structure in the model:”
plot (mixto1pc1, Dim.1 ~ fitted(.) | Postura, abline = c (0,1))
Efecto aleatorio: La caja
plot(mixto1pc1, resid(., scaled=TRUE) ~ fitted(.) |Caja,abline = 0)
#The following output shows the scaled residuals for each group:
plot(mixto1pc1, Caja ~ resid(., scaled=TRUE)) 
plot (mixto1pc1, Dim.1 ~ fitted(.) | Caja, abline = c (0,1))
3. Calculando la magnitud del efecto del modelo
La magnitud del efecto:
Effect sizes, are metrics that represent the amount of differences between two sample means.
Cohen’s d and Hedges’ g are interpreted in a similar way. Cohen suggested using the following rule of thumb for interpreting results:
Small effect (cannot be discerned by the naked eye) = 0.2 Medium Effect = 0.5 Large Effect (can be seen by the naked eye) = 0.8
Cohen, J. (1977). Statistical power analysis for the behavioral sciences. Routledge.
library(effectsize)
#proporcion de varianza atribuida a una variable
eta_squared(mixto1pc1)## # Effect Size for ANOVA (Type III)
##
## Parameter | Eta2 (partial) | 95% CI
## -----------------------------------------
## TREATMENT | 0.47 | [0.07, 1.00]
##
## - One-sided CIs: upper bound fixed at [1.00].cohens_f(mixto1pc1)## # Effect Size for ANOVA (Type III)
##
## Parameter | Cohen's f (partial) | 95% CI
## ---------------------------------------------
## TREATMENT | 0.95 | [0.27, Inf]
##
## - One-sided CIs: upper bound fixed at [Inf].Para contrastar entre tratamientos:
library(emmeans)
trat.emm.s <- emmeans(mixto1pc1, "TREATMENT")
pairs(trat.emm.s)## contrast estimate SE df t.ratio p.value
## DD - HW 0.594 0.258 11.7 2.306 0.0943
## DD - LW 0.834 0.260 12.1 3.210 0.0188
## HW - LW 0.240 0.260 12.2 0.924 0.6362
##
## Degrees-of-freedom method: kenward-roger
## P value adjustment: tukey method for comparing a family of 3 estimatesA partir de esos contrastes, la magnitud del efecto:
eff_size(trat.emm.s, sigma = sigma(mixto1pc1), edf = 12)## contrast effect.size SE df lower.CL upper.CL
## DD - HW 0.820 0.367 11.3 0.0157 1.62
## DD - LW 1.152 0.359 11.3 0.3636 1.94
## HW - LW 0.332 0.365 11.4 -0.4692 1.13
##
## sigma used for effect sizes: 0.7245
## Degrees-of-freedom method: inherited from kenward-roger when re-gridding
## Confidence level used: 0.95Plotting Random Effects of Mixed Models
Visualizing (generalized) linear mixed effects models
Comparación con DD
library(sjPlot)
library(lme4)
datos_morfo_trans <- read.table("morfo_med_corr.txt",h=T)
#para volver una variable categorica:
datos_morfo_trans$TREATMENT <- as.factor(datos_morfo_trans$TREATMENT)
datos_morfo_trans$TREATMENT <-relevel(datos_morfo_trans$TREATMENT, ref="DD")
mixto1pc1_DD<-lmer(Pcrot~TREATMENT+(1|Caja)+(1|Postura),data=datos_morfo_trans)
plot_model(mixto1pc1_DD,colors = c("#33638DFF"),vline.color ="black")
Comparación con HW
#para volver una variable categorica:
datos_morfo_trans$TREATMENT <- as.factor(datos_morfo_trans$TREATMENT)
datos_morfo_trans$TREATMENT<-relevel(datos_morfo_trans$TREATMENT, ref="HW")
mixto1pc1_HW<-lmer(Pcrot~TREATMENT+(1|Caja)+(1|Postura),data=datos_morfo_trans)
summary(mixto1pc1_HW)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Pcrot ~ TREATMENT + (1 | Caja) + (1 | Postura)
## Data: datos_morfo_trans
##
## REML criterion at convergence: 421.9
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.57051 -0.62221 -0.03631 0.61635 2.54916
##
## Random effects:
## Groups Name Variance Std.Dev.
## Caja (Intercept) 0.1713 0.4139
## Postura (Intercept) 0.2663 0.5161
## Residual 0.5249 0.7245
## Number of obs: 175, groups: Caja, 21; Postura, 7
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -0.1370 0.2670 11.3731 -0.513 0.6176
## TREATMENTDD 0.5940 0.2576 11.9251 2.306 0.0399 *
## TREATMENTLW -0.2403 0.2601 12.3826 -0.924 0.3732
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) TREATMENTD
## TREATMENTDD -0.483
## TREATMENTLW -0.479 0.496plot_model(mixto1pc1_HW,colors = c("#33638DFF"))
Comparación con LW
#para volver una variable categorica:
datos_morfo_trans$TREATMENT <- as.factor(datos_morfo_trans$TREATMENT)
datos_morfo_trans$TREATMENT <-relevel(datos_morfo_trans$TREATMENT, ref="LW")
mixto1pc1_LW<-lmer(Pcrot~TREATMENT+(1|Caja)+(1|Postura),data=datos_morfo_trans)
summary(mixto1pc1_LW)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Pcrot ~ TREATMENT + (1 | Caja) + (1 | Postura)
## Data: datos_morfo_trans
##
## REML criterion at convergence: 421.9
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.57051 -0.62221 -0.03631 0.61635 2.54916
##
## Random effects:
## Groups Name Variance Std.Dev.
## Caja (Intercept) 0.1713 0.4139
## Postura (Intercept) 0.2663 0.5161
## Residual 0.5249 0.7245
## Number of obs: 175, groups: Caja, 21; Postura, 7
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -0.3774 0.2692 11.7403 -1.402 0.18684
## TREATMENTHW 0.2403 0.2601 12.3826 0.924 0.37318
## TREATMENTDD 0.8344 0.2598 12.3280 3.211 0.00725 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) TREATMENTH
## TREATMENTHW -0.492
## TREATMENTDD -0.492 0.509plot_model(mixto1pc1_LW,colors = c("#33638DFF"))
Probando otras opciones para ver magnitud del efecto:
library(effectsize)
effectsize(mixto1pc1_DD)## # Standardization method: refit
##
## Parameter | Std. Coef. | 95% CI
## -----------------------------------------
## (Intercept) | 0.46 | [-0.07, 0.98]
## TREATMENTHW | -0.59 | [-1.10, -0.09]
## TREATMENTLW | -0.83 | [-1.35, -0.32]effectsize(mixto1pc1_HW)## # Standardization method: refit
##
## Parameter | Std. Coef. | 95% CI
## ----------------------------------------
## (Intercept) | -0.14 | [-0.66, 0.39]
## TREATMENTDD | 0.59 | [ 0.09, 1.10]
## TREATMENTLW | -0.24 | [-0.75, 0.27]effectsize(mixto1pc1_LW)## # Standardization method: refit
##
## Parameter | Std. Coef. | 95% CI
## ----------------------------------------
## (Intercept) | -0.38 | [-0.91, 0.15]
## TREATMENTHW | 0.24 | [-0.27, 0.75]
## TREATMENTDD | 0.83 | [ 0.32, 1.35]El modelo mixto con PC2
datos_morfo_trans$TREATMENT <- as.factor(datos_morfo_trans$TREATMENT)
datos_morfo_trans$TREATMENT <-relevel(datos_morfo_trans$TREATMENT, ref="DD")
mixto1pc2<-lmer(Dim.2~TREATMENT+(1|Caja)+(1|Postura),data=datos_morfo_trans)
summary(mixto1pc2)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Dim.2 ~ TREATMENT + (1 | Caja) + (1 | Postura)
## Data: datos_morfo_trans
##
## REML criterion at convergence: 301.9
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.33675 -0.67522 -0.03272 0.60761 2.69881
##
## Random effects:
## Groups Name Variance Std.Dev.
## Caja (Intercept) 0.01687 0.1299
## Postura (Intercept) 0.09440 0.3072
## Residual 0.28331 0.5323
## Number of obs: 175, groups: Caja, 21; Postura, 7
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -0.06993 0.14334 9.48768 -0.488 0.637
## TREATMENTLW 0.10352 0.12161 12.58436 0.851 0.411
## TREATMENTHW 0.10150 0.11917 11.70540 0.852 0.411
##
## Correlation of Fixed Effects:
## (Intr) TREATMENTL
## TREATMENTLW -0.405
## TREATMENTHW -0.413 0.487library(sjPlot)
library(lme4)
tab_model(mixto1pc2)| Dim.2 | |||
|---|---|---|---|
| Predictors | Estimates | CI | p |
| (Intercept) | -0.07 | -0.35 – 0.21 | 0.626 |
| TREATMENT [LW] | 0.10 | -0.14 – 0.34 | 0.396 |
| TREATMENT [HW] | 0.10 | -0.13 – 0.34 | 0.396 |
| Random Effects | |||
| σ2 | 0.28 | ||
| τ00 Caja | 0.02 | ||
| τ00 Postura | 0.09 | ||
| ICC | 0.28 | ||
| N Caja | 21 | ||
| N Postura | 7 | ||
| Observations | 175 | ||
| Marginal R2 / Conditional R2 | 0.006 / 0.286 | ||
anova(mixto1pc2)## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## TREATMENT 0.27623 0.13812 2 12.301 0.4875 0.6255Haciendo el ajuste en la significancia:
tab_model(mixto1pc2, p.val = "kr", show.df = TRUE)| Dim.2 | ||||
|---|---|---|---|---|
| Predictors | Estimates | CI | p | df |
| (Intercept) | -0.07 | -0.39 – 0.25 | 0.637 | 9.55 |
| TREATMENT [LW] | 0.10 | -0.16 – 0.37 | 0.411 | 12.21 |
| TREATMENT [HW] | 0.10 | -0.16 – 0.36 | 0.412 | 11.35 |
| Random Effects | ||||
| σ2 | 0.28 | |||
| τ00 Caja | 0.02 | |||
| τ00 Postura | 0.09 | |||
| ICC | 0.28 | |||
| N Caja | 21 | |||
| N Postura | 7 | |||
| Observations | 175 | |||
| Marginal R2 / Conditional R2 | 0.006 / 0.286 | |||
##RESULTADOS DESEMPEÑO LOCOMOTOR
Para este experimento se utilizaron 7 posturas diferentes. Cada postura se dividió en 3 tratamientos (Agua Baja AB,Agua Estable AE ,Agua Desecación AD), para un total de 21 cajas.5 individuos provenientes de cada caja fueron utilizados para las pruebas de desempeño locomotor. Estas pruebas se realizaron entre las 9 y las 16 horas del día. Se evaluó el desempeño locomotor a 25, 30, 33, 37 y 40 grados centígrados. El desempeño de cada renacuajo fue evaluado a una única temperatura y para esto el animal fue introducido en recipiente con agua durante 3 minutos y se filmó su respuesta al ser estimulado con un estilete. A partir de los videos y utlizando el programa Tracker se determinó para cada renacuajo los valores máximos de velocidad, aceleración y distancia.
programa para análisis de videos
knitr::include_graphics("locom.png")Los datos:
desloc<- read.table("desloc_def.txt",h=T)
head(desloc)## TREATMENT postura caja temp Id_Random_Foto individuo vel acel
## 1 HW C2 16 37.0 1 4 0.2733747 3.651805
## 2 HW C3 17 25.5 2 1 0.2534460 2.922191
## 3 DD C5 12 33.0 3 5 0.3889761 4.750895
## 4 HW C2 16 40.0 4 6 0.1649230 1.672477
## 5 DD C1 8 40.0 5 6 0.3145094 2.745355
## 6 DD C6 13 30.0 6 3 0.5120589 6.788423
## long vel_a mag_p peso Longitud_Total Peso Altura_Cabeza
## 1 0.09981042 121.49333 0.10714760 0.030 24.97633 0.030 2.486333
## 2 0.08462006 167.03079 0.09233168 0.042 23.77367 0.042 2.164333
## 3 0.80143725 288.22802 0.17048985 0.035 26.36967 0.035 2.817333
## 4 0.16512646 53.92223 0.12073766 0.017 18.02033 0.017 1.881667
## 5 0.06470735 554.99910 0.08438124 0.070 29.31800 0.070 2.485333
## 6 0.39807844 186.40142 0.16396688 0.021 20.78933 0.021 2.109000
## Altura_Cola Longitud_Cola Perim_Musculo1 Longitud_Cabeza L_Longitud_Total
## 1 1.565333 17.52700 38.31700 7.633000 1.397529
## 2 1.467667 16.84500 34.71633 7.114333 1.376096
## 3 1.675667 18.74700 38.77133 7.669667 1.421105
## 4 1.168667 12.53900 27.13233 5.639667 1.255763
## 5 2.054333 20.18367 43.89933 8.744000 1.467134
## 6 1.413333 15.22967 32.12367 5.471000 1.317841
## L_Peso L_Altura_Cabeza L_Altura_Cola L_Longitud_Cola L_Perimetro_Musculo
## 1 -1.522879 0.3955594 0.19460683 1.243708 1.583391
## 2 -1.376751 0.3353241 0.16662743 1.226471 1.540534
## 3 -1.455932 0.4498382 0.22418763 1.272932 1.588511
## 4 -1.769551 0.2745427 0.06769066 1.098263 1.433487
## 5 -1.154902 0.3953846 0.31267091 1.305000 1.642458
## 6 -1.677781 0.3240766 0.15024460 1.182690 1.506825
## L_Longitud_Cabeza Dim.1 Dim.2 Pcrot
## 1 0.8826953 3.2632 0.2259 1.4708136
## 2 0.8521342 2.4854 -0.4660 1.1202452
## 3 0.8847765 4.1504 0.5298 1.8707077
## 4 0.7512534 -1.4971 0.3299 -0.6747671
## 5 0.9417101 5.8749 -0.5259 2.6480109
## 6 0.7380667 0.4226 -0.1818 0.1904794Librerias para las figuras:
library(ggplot2)
library(yarrr)
library(ggstance)
library(ggformula)
library(viridisLite)
library(viridis)Figuras exploratorias
#####Aceleración
desloc<- read.table("desloc_def.txt",h=T)
desloc$TREATMENT <- as.factor(desloc$TREATMENT)
desloc$TREATMENT <-relevel(desloc$TREATMENT, ref="DD")ggplot(desloc,aes(x=temp, y=acel,color=TREATMENT)) +
geom_point() + geom_smooth(aes(fill=TREATMENT),span=5 ) + scale_color_manual(values = c( "#33638DFF","#B8DE29FF","#FDE725FF"))+ theme_bw() +scale_fill_manual(values = c("#33638DFF","#B8DE29FF","#FDE725FF"))+labs(x= "Temperature (°C)", y = "Maximum acceleration (m/s^2)")+theme(plot.title = element_text(size = 14, family = "Tahoma", face = "bold"), text = element_text(size = 14, family ="Tahoma"),axis.title = element_text(face="bold"),axis.text.x=element_text(size = 12))
###Modelo MixtoGAM
Paquetes para hacer el análisis
library(lme4)
library(mgcv)
library(gamm4)
library(sjPlot)desloc$TREATMENT <-relevel(desloc$TREATMENT, ref="DD")
mod1.gam = gam(acel ~ s(temp,k=5)+TREATMENT, data=desloc)
summary(mod1.gam)##
## Family: gaussian
## Link function: identity
##
## Formula:
## acel ~ s(temp, k = 5) + TREATMENT
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.1384 0.1668 24.814 < 2e-16 ***
## TREATMENTHW -0.8273 0.2358 -3.508 0.000682 ***
## TREATMENTLW -0.7143 0.2376 -3.006 0.003357 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(temp) 2.35 2.83 7.54 0.000196 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.251 Deviance explained = 28.2%
## GCV = 1.0262 Scale est. = 0.9734 n = 104anova(mod1.gam)##
## Family: gaussian
## Link function: identity
##
## Formula:
## acel ~ s(temp, k = 5) + TREATMENT
##
## Parametric Terms:
## df F p-value
## TREATMENT 2 7.214 0.00119
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(temp) 2.35 2.83 7.54 0.000196plot_model(mod1.gam,colors = c("#33638DFF"))
tab_model(mod1.gam, p.val = "kr", show.df = TRUE)| acel | ||||
|---|---|---|---|---|
| Predictors | Estimates | CI | p | df |
| (Intercept) | 4.14 | 3.81 – 4.47 | <0.001 | 98.65 |
| TREATMENT [HW] | -0.83 | -1.30 – -0.36 | 0.001 | 98.65 |
| TREATMENT [LW] | -0.71 | -1.19 – -0.24 | 0.003 | 98.65 |
| Smooth term (temp) | <0.001 | |||
| Observations | 104 | |||
| R2 | 0.251 | |||
Para contrastar entre tratamientos:
library(emmeans)
trat.emm.s <- emmeans(mod1.gam, "TREATMENT")
pairs(trat.emm.s)## contrast estimate SE df t.ratio p.value
## DD - HW 0.827 0.236 98.7 3.508 0.0020
## DD - LW 0.714 0.238 98.7 3.006 0.0093
## HW - LW -0.113 0.238 98.7 -0.475 0.8830
##
## P value adjustment: tukey method for comparing a family of 3 estimatesdesloc$tTREATMENT <- as.factor(desloc$TREATMENT)
desloc$TREATMENT <-relevel(desloc$TREATMENT, ref="HW")
mod1.gam_HW = gam(acel ~ s(temp,k=5)+TREATMENT, data=desloc)
summary(mod1.gam_HW)##
## Family: gaussian
## Link function: identity
##
## Formula:
## acel ~ s(temp, k = 5) + TREATMENT
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.3112 0.1668 19.854 < 2e-16 ***
## TREATMENTDD 0.8273 0.2358 3.508 0.000682 ***
## TREATMENTLW 0.1130 0.2376 0.475 0.635511
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(temp) 2.35 2.83 7.54 0.000196 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.251 Deviance explained = 28.2%
## GCV = 1.0262 Scale est. = 0.9734 n = 104anova(mod1.gam_HW)##
## Family: gaussian
## Link function: identity
##
## Formula:
## acel ~ s(temp, k = 5) + TREATMENT
##
## Parametric Terms:
## df F p-value
## TREATMENT 2 7.214 0.00119
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(temp) 2.35 2.83 7.54 0.000196plot_model(mod1.gam_HW)
tab_model(mod1.gam_HW, p.val = "kr", show.df = TRUE)| acel | ||||
|---|---|---|---|---|
| Predictors | Estimates | CI | p | df |
| (Intercept) | 3.31 | 2.98 – 3.64 | <0.001 | 98.65 |
| TREATMENT [DD] | 0.83 | 0.36 – 1.30 | 0.001 | 98.65 |
| TREATMENT [LW] | 0.11 | -0.36 – 0.58 | 0.636 | 98.65 |
| Smooth term (temp) | <0.001 | |||
| Observations | 104 | |||
| R2 | 0.251 | |||
library(effectsize)
effectsize(mod1.gam)## # Standardization method: refit
##
## Component | Parameter | Std. Coef. | 95% CI
## --------------------------------------------------------
## conditional | (Intercept) | 0.45 | [ 0.16, 0.74]
## conditional | TREATMENTHW | -0.73 | [-1.14, -0.32]
## conditional | TREATMENTLW | -0.63 | [-1.04, -0.21]
## smooth_terms | s(temp) | |desloc$TREATMENT <-relevel(desloc$TREATMENT, ref="DD")
mod2.gam = gam(acel ~ s(temp,k=5)+TREATMENT+postura, data=desloc)
summary(mod2.gam)##
## Family: gaussian
## Link function: identity
##
## Formula:
## acel ~ s(temp, k = 5) + TREATMENT + postura
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.23104 0.28073 15.072 < 2e-16 ***
## TREATMENTHW -0.82729 0.22906 -3.612 0.000494 ***
## TREATMENTLW -0.71910 0.23091 -3.114 0.002454 **
## posturaC2 0.16988 0.34990 0.486 0.628462
## posturaC3 -0.24780 0.35644 -0.695 0.488661
## posturaC4 -0.81810 0.34990 -2.338 0.021534 *
## posturaC5 0.15368 0.34990 0.439 0.661531
## posturaC6 0.07425 0.34990 0.212 0.832403
## posturaC7 0.02050 0.34990 0.059 0.953403
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(temp) 2.373 2.855 8.054 0.000113 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.293 Deviance explained = 36.4%
## GCV = 1.031 Scale est. = 0.91822 n = 104anova(mod2.gam)##
## Family: gaussian
## Link function: identity
##
## Formula:
## acel ~ s(temp, k = 5) + TREATMENT + postura
##
## Parametric Terms:
## df F p-value
## TREATMENT 2 7.682 0.000819
## postura 6 1.981 0.076269
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(temp) 2.373 2.855 8.054 0.000113plot_model(mod2.gam)
tab_model(mod2.gam, p.val = "kr", show.df = TRUE)| acel | ||||
|---|---|---|---|---|
| Predictors | Estimates | CI | p | df |
| (Intercept) | 4.23 | 3.67 – 4.79 | <0.001 | 92.63 |
| TREATMENT [HW] | -0.83 | -1.28 – -0.37 | <0.001 | 92.63 |
| TREATMENT [LW] | -0.72 | -1.18 – -0.26 | 0.002 | 92.63 |
| postura [C2] | 0.17 | -0.52 – 0.86 | 0.628 | 92.63 |
| postura [C3] | -0.25 | -0.96 – 0.46 | 0.489 | 92.63 |
| postura [C4] | -0.82 | -1.51 – -0.12 | 0.022 | 92.63 |
| postura [C5] | 0.15 | -0.54 – 0.85 | 0.662 | 92.63 |
| postura [C6] | 0.07 | -0.62 – 0.77 | 0.832 | 92.63 |
| postura [C7] | 0.02 | -0.67 – 0.72 | 0.953 | 92.63 |
| Smooth term (temp) | <0.001 | |||
| Observations | 104 | |||
| R2 | 0.293 | |||
Desempeño locomotor con el tamaño corporal como covariable
El componente principal rotado correspondiente a Body size (Pcrot) se usó como covariable del modelo
desloc$TREATMENT <-relevel(desloc$TREATMENT, ref="DD")
mod3.gam = gam(acel ~ s(temp,k=5)+TREATMENT+Pcrot, data=desloc)
summary(mod3.gam)##
## Family: gaussian
## Link function: identity
##
## Formula:
## acel ~ s(temp, k = 5) + TREATMENT + Pcrot
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.93881 0.17346 22.707 < 2e-16 ***
## TREATMENTHW -0.61372 0.23749 -2.584 0.01124 *
## TREATMENTLW -0.50831 0.23841 -2.132 0.03551 *
## Pcrot 0.29687 0.09863 3.010 0.00332 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(temp) 2.411 2.897 7.633 0.000149 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.308 Deviance explained = 34.4%
## GCV = 0.95797 Scale est. = 0.89892 n = 104anova(mod3.gam)##
## Family: gaussian
## Link function: identity
##
## Formula:
## acel ~ s(temp, k = 5) + TREATMENT + Pcrot
##
## Parametric Terms:
## df F p-value
## TREATMENT 2 3.723 0.02765
## Pcrot 1 9.060 0.00332
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(temp) 2.411 2.897 7.633 0.000149plot_model(mod3.gam,colors = c("#33638DFF"))
tab_model(mod3.gam, p.val = "kr", show.df = TRUE)| acel | ||||
|---|---|---|---|---|
| Predictors | Estimates | CI | p | df |
| (Intercept) | 3.94 | 3.59 – 4.28 | <0.001 | 97.59 |
| TREATMENT [HW] | -0.61 | -1.09 – -0.14 | 0.011 | 97.59 |
| TREATMENT [LW] | -0.51 | -0.98 – -0.04 | 0.036 | 97.59 |
| Pcrot | 0.30 | 0.10 – 0.49 | 0.003 | 97.59 |
| Smooth term (temp) | <0.001 | |||
| Observations | 104 | |||
| R2 | 0.308 | |||
desloc$TREATMENT <-relevel(desloc$TREATMENT, ref="HW")
mod4.gam = gam(acel ~ s(temp,k=5)+TREATMENT*Pcrot, data=desloc)
summary(mod4.gam)##
## Family: gaussian
## Link function: identity
##
## Formula:
## acel ~ s(temp, k = 5) + TREATMENT * Pcrot
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.33567 0.15951 20.911 < 2e-16 ***
## TREATMENTDD 0.71846 0.24411 2.943 0.00408 **
## TREATMENTLW 0.09687 0.22713 0.427 0.67069
## Pcrot 0.51757 0.18993 2.725 0.00765 **
## TREATMENTDD:Pcrot -0.39247 0.23548 -1.667 0.09885 .
## TREATMENTLW:Pcrot -0.11569 0.27522 -0.420 0.67516
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(temp) 2.43 2.917 7.708 0.000128 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.317 Deviance explained = 36.6%
## GCV = 0.96593 Scale est. = 0.88763 n = 104anova(mod4.gam)##
## Family: gaussian
## Link function: identity
##
## Formula:
## acel ~ s(temp, k = 5) + TREATMENT * Pcrot
##
## Parametric Terms:
## df F p-value
## TREATMENT 2 4.870 0.00968
## Pcrot 1 7.426 0.00765
## TREATMENT:Pcrot 2 1.581 0.21114
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(temp) 2.430 2.917 7.708 0.000128plot_model(mod4.gam,colors = c("#33638DFF"))
tab_model(mod4.gam, p.val = "kr", show.df = TRUE)| acel | ||||
|---|---|---|---|---|
| Predictors | Estimates | CI | p | df |
| (Intercept) | 3.34 | 3.02 – 3.65 | <0.001 | 95.57 |
| TREATMENT [DD] | 0.72 | 0.23 – 1.20 | 0.004 | 95.57 |
| TREATMENT [LW] | 0.10 | -0.35 – 0.55 | 0.671 | 95.57 |
| Pcrot | 0.52 | 0.14 – 0.89 | 0.008 | 95.57 |
| TREATMENT [DD] × Pcrot | -0.39 | -0.86 – 0.07 | 0.099 | 95.57 |
| TREATMENT [LW] × Pcrot | -0.12 | -0.66 – 0.43 | 0.675 | 95.57 |
| Smooth term (temp) | <0.001 | |||
| Observations | 104 | |||
| R2 | 0.317 | |||
library(effectsize)
effectsize(mod3.gam)## # Standardization method: refit
##
## Component | Parameter | Std. Coef. | 95% CI
## --------------------------------------------------------
## conditional | (Intercept) | 0.33 | [ 0.04, 0.62]
## conditional | TREATMENTHW | -0.54 | [-0.95, -0.12]
## conditional | TREATMENTLW | -0.45 | [-0.86, -0.03]
## conditional | Pcrot | 0.26 | [ 0.09, 0.44]
## smooth_terms | s(temp) | |Errores modelo mixto-GAM
library(mgcv) library(lme4) desloc1\(tratamiento <- as.factor(desloc1\)tratamiento) desloc1\(postura <- as.factor(desloc1\)postura) mixtogam<-gamm(acel~tratamiento+s(temp,k=3)+(1|postura),data=desloc1)
Error in gamm(acel ~ tratamiento + s(temp, k = 3) + (1 | postura), data = desloc1) : Not enough (non-NA) data to do anything meaningful
mixtogamo<-gamm(acel(TREATMENT)+s(temp),random=list(postura=1), data=desloc) Error in smooth.construct.tp.smooth.spec(object, dk\(data, dk\)knots) : A term has fewer unique covariate combinations than specified maximum degrees of freedom