Absorción de dióxido de carbono en plantas de pasto:
Se midió la absorción de CO2 de la especie de pasto Echinochloa crus-galliseis, (seis plantas de Quebec y seis plantas de Mississippi) a varios niveles de concentración ambiental(95,175,250,350,500,675,1000) (mL/L). La mitad de las plantas de cada tipo se enfriaron durante la noche antes de realizar el experimento y la otra mitad no(nonchilled-chilled). Se midió la absorción (μmol/m2) de cada en planta con tratamiento nonchilled y chilled en cada nivel de concentración.
library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.2 ✔ readr 2.1.4
## ✔ forcats 1.0.0 ✔ stringr 1.5.0
## ✔ ggplot2 3.4.2 ✔ tibble 3.2.1
## ✔ lubridate 1.9.2 ✔ tidyr 1.3.0
## ✔ purrr 1.0.1
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(broom.mixed)
library(lme4)
## Loading required package: Matrix
##
## Attaching package: 'Matrix'
##
## The following objects are masked from 'package:tidyr':
##
## expand, pack, unpack
data("CO2")
uptake2 <- CO2$uptake + runif(84, 0.02,0.04)
#View(CO2)
ggplot(CO2, aes(x = conc, y = uptake2, color = Type)) + geom_point(aes(shape = Treatment)) +
geom_path(aes(group = Plant, lty = Treatment)) +
theme_bw()
Se genera un gráfico de dispersión donde los puntos están coloreados
según el tipo de planta (Type) y las formas de los puntos representan
los diferentes tratamientos (Treatment). Además, se agregan líneas
conectando los puntos de cada planta (Plant) según el tratamiento.
La gráfica muestra que la mayor absorción de CO2 se obtiene de las plantas de Quebec, con y sin tratamiento de enfriamiento. Las plantas de Mississippi siempre se mantienen abajo en en cuanto al consumo de CO2
En este ejercicio el efecto aleatorio es la planta (Plant, debido a que captura variaciones entre cada unidad observada, estas variaciones no se explican por los efectos fijos, Tipo (Type), Tratamiento (Trantment) o la cocentración de CO2 (conc)
#la grafica muestra que hay interacción entre los tratamiento
Fit1 <- lm(uptake2 ~ I(log(conc)) + Type:Treatment, data = CO2) #modelo lineal simple
Fit2 <- lmer(uptake2 ~ I(log(conc)) + Type:Treatment + (1 | Plant), data = CO2) # modelo lineal mixto (se establece como factor aleatorio a la Planta (Plant))
## fixed-effect model matrix is rank deficient so dropping 1 column / coefficient
broom::glance(Fit1)
## # A tibble: 1 × 12
## r.squared adj.r.squared sigma statistic p.value df logLik AIC BIC
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 0.823 0.814 4.67 91.6 7.23e-29 4 -246. 504. 519.
## # ℹ 3 more variables: deviance <dbl>, df.residual <int>, nobs <int>
broom::tidy(Fit1)
## # A tibble: 6 × 5
## term estimate std.error statistic p.value
## <chr> <dbl> <dbl> <dbl> <dbl>
## 1 (Intercept) -33.5 4.08 -8.22 3.25e-12
## 2 I(log(conc)) 8.48 0.679 12.5 2.07e-20
## 3 TypeQuebec:Treatmentnonchilled 19.5 1.44 13.5 2.69e-22
## 4 TypeMississippi:Treatmentnonchilled 10.1 1.44 7.03 6.42e-10
## 5 TypeQuebec:Treatmentchilled 15.9 1.44 11.1 1.02e-17
## 6 TypeMississippi:Treatmentchilled NA NA NA NA
broom.mixed::glance(Fit2)
## # A tibble: 1 × 7
## nobs sigma logLik AIC BIC REMLcrit df.residual
## <int> <dbl> <dbl> <dbl> <dbl> <dbl> <int>
## 1 84 4.50 -241. 496. 513. 482. 77
broom.mixed::tidy(Fit2)
## # A tibble: 7 × 6
## effect group term estimate std.error statistic
## <chr> <chr> <chr> <dbl> <dbl> <dbl>
## 1 fixed <NA> (Intercept) -33.5 4.02 -8.33
## 2 fixed <NA> I(log(conc)) 8.48 0.654 13.0
## 3 fixed <NA> TypeQuebec:Treatmentnonchilled 19.5 1.84 10.6
## 4 fixed <NA> TypeMississippi:Treatmentnonch… 10.1 1.84 5.52
## 5 fixed <NA> TypeQuebec:Treatmentchilled 15.9 1.84 8.68
## 6 ran_pars Plant sd__(Intercept) 1.47 NA NA
## 7 ran_pars Residual sd__Observation 4.50 NA NA
Se ajustan dos modelos estadísticos, un modelo lineal (Fit1) y un modelo lineal mixto (Fit2), se utilizan las funciones glance y tidy para obtener resúmenes del ajuste de los modelos Fit1 y Fit2, respectivamente. Broom.mixed se utiliza para obtener un resumen específico del modelo lineal mixto Fit2.
Como se ha establecido a la var
broom.mixed::tidy(Fit2) %>% head() #estimador para cada planta y para cada observacion
## # A tibble: 6 × 6
## effect group term estimate std.error statistic
## <chr> <chr> <chr> <dbl> <dbl> <dbl>
## 1 fixed <NA> (Intercept) -33.5 4.02 -8.33
## 2 fixed <NA> I(log(conc)) 8.48 0.654 13.0
## 3 fixed <NA> TypeQuebec:Treatmentnonchilled 19.5 1.84 10.6
## 4 fixed <NA> TypeMississippi:Treatmentnonchill… 10.1 1.84 5.52
## 5 fixed <NA> TypeQuebec:Treatmentchilled 15.9 1.84 8.68
## 6 ran_pars Plant sd__(Intercept) 1.47 NA NA
Como se ha establecido a la variable planta (Plant) como efecto aleatorio, nos fijamos en los efectos fijos para confirmar qué interacción tiene el mejor resultado en cuanto a la medición de absorción de CO2.
Conclución: Como se había observado en la gráfica, la mejor interacción la ofrece TypeQuebec:Treatmentnonchilled, es decir, plantas de Quebec las cuales se enfriaron antes de realizar el experimento, el estimado de esta observación es de 19.513360. La segunda mejor estimación la ofrece TypeQuebec:Treatmentchilled con una valor de 15.930936, esto confirma que las plantas de Quebec tienen una mayor absorción de CO2, independientemente de que tratamiento se usó antes de medir la absorción.