knitr::opts_chunk$set(echo = TRUE)
library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.4 ✔ readr 2.1.5
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ ggplot2 3.5.1 ✔ tibble 3.2.1
## ✔ lubridate 1.9.3 ✔ tidyr 1.3.1
## ✔ purrr 1.0.2
## ── 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(knitr)
library(data.table)
## Warning: package 'data.table' was built under R version 4.4.3
##
## Adjuntando el paquete: 'data.table'
##
## The following objects are masked from 'package:lubridate':
##
## hour, isoweek, mday, minute, month, quarter, second, wday, week,
## yday, year
##
## The following objects are masked from 'package:dplyr':
##
## between, first, last
##
## The following object is masked from 'package:purrr':
##
## transpose
library(haven)
library(geepack)
## Warning: package 'geepack' was built under R version 4.4.3
library(corrplot)
## Warning: package 'corrplot' was built under R version 4.4.3
## corrplot 0.95 loaded
library(psych)
##
## Adjuntando el paquete: 'psych'
##
## The following objects are masked from 'package:ggplot2':
##
## %+%, alpha
#library(sjPlot) opcional para ver resultados = tab_model(modelo)
library(kableExtra)
## Warning: package 'kableExtra' was built under R version 4.4.3
##
## Adjuntando el paquete: 'kableExtra'
##
## The following object is masked from 'package:dplyr':
##
## group_rows
options(scipen=999, digits=3)
CD4 <- read_dta("C:/Users/Administrador/Downloads/CD4.dta")
#View(CD4_Ejercicio1)
colSums(is.na(CD4)) ### Valores Faltantes
## id time cd4 Age gender sex Treatment
## 0 0 42 0 0 0 0
La base esta en formato “LONG” Los datos se muestran en formato “long”, donde cada fila es una medición y time es el momento de la medición. El id se repite para cada individuo de acuerdo al n de mediciones que tuvo. 1)Evaluar el cambio de CD4 en el tiempo ajustando por edad y sexo como posibles confundidores 2) Evaluar el cambio de CD4 en función del tiempo EN HOMBRES Y MUJERES 3) Evaluar el cambio de CD4 en función del tiempo en los grupos de TRATAMIENTO
CD4$Treatment <- as.numeric(CD4$Treatment)
# LINK: https://tidyr.tidyverse.org/articles/pivot.html#wider
# Para ver las correlaciones, debemos pasar lods datos a formato "WIDE"
CD4_WIDE <- CD4 %>%
pivot_wider(names_from = time, ## selecciono la variable repetida
names_prefix = "Tiempo", ## agrego el prefijo para el nuevo nombre
values_from = cd4) ### selecciono de donde va a tomar los valores
# AL PASAR DE LONG A WIDE se ELIMINA LA VARIABLE TIME Y SE CREAN LAS 3
# MEDIDAS REPETIDAS DE CD4 CON LA INFORMACION DE TIME, Tiempo0, Tiempo6 Y Tiempo12
CD4_WIDE[1:10,-3]
## # A tibble: 10 × 7
## id Age sex Treatment Tiempo12 Tiempo6 Tiempo0
## <dbl+lbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 18 [DIC-107] 23 1 0 529 281 364
## 2 47 [JB-11] 30 1 0 506 605 311
## 3 11 [BKC-32] 33 1 0 735 780 554
## 4 71 [MLL-61*] 47 1 0 670 638 468
## 5 69 [ML-37] 44 1 0 365 487 281
## 6 65 [MDL-81] 38 0 0 1024 850 704
## 7 14 [CCG-45] 35 1 0 613 524 377
## 8 75 [NEC-119] 31 1 0 729 590 521
## 9 20 [DK-04] 34 1 0 682 777 652
## 10 62 [LFD-39] 26 1 0 434 372 607
View(CD4_WIDE)
head(CD4[order(CD4$id, decreasing = F),]) %>% kable()
| id | time | cd4 | Age | gender | sex | Treatment |
|---|---|---|---|---|---|---|
| 1 | 0 | 311 | 62 | M | 1 | 1 |
| 1 | 6 | 473 | 62 | M | 1 | 1 |
| 1 | 12 | 408 | 62 | M | 1 | 1 |
| 2 | 0 | 436 | 25 | M | 1 | 0 |
| 2 | 6 | 442 | 25 | M | 1 | 0 |
| 2 | 12 | 643 | 25 | M | 1 | 0 |
view(CD4)
corrplot.mixed(cor(CD4_WIDE[,6:8], # selecciono las variables dependientes correlacionadas
use = "pairwise.complete"), # correlacion
lower = "number", # panel inferior = número
upper = "circle", # panel superior = Círculos
title = "Correlacion CD4") # título
CD4_Long <- CD4_WIDE %>% # selecciono base wide
pivot_longer(
cols = starts_with("Tiempo"), # columnas que empiezan con tiempo (0, 6 ,12)
names_to = "tiempo_nuevo", # Nueba variable que va a contener tiempo(0, 6, 12)
names_transform = as.integer, # transformo tiempo a integer-numerica
names_prefix = "Tiempo", # saco la palabra tiempo para que no quede tiempo6 como valor
values_to = c("cd4_nuevo") # los valores de cd4 de cada tiempo a la nueva variable
)
head(CD4_Long)
## # A tibble: 6 × 7
## id Age gender sex Treatment tiempo_nuevo cd4_nuevo
## <dbl+lbl> <dbl> <chr> <dbl> <dbl> <int> <dbl>
## 1 18 [DIC-107] 23 M 1 0 12 529
## 2 18 [DIC-107] 23 M 1 0 6 281
## 3 18 [DIC-107] 23 M 1 0 0 364
## 4 47 [JB-11] 30 M 1 0 12 506
## 5 47 [JB-11] 30 M 1 0 6 605
## 6 47 [JB-11] 30 M 1 0 0 311
#modelo de regresión lineal si las variables fueran independientes
mod_lineal <- lm(cd4~time +Age + sex, data = CD4)
summary(mod_lineal)
##
## Call:
## lm(formula = cd4 ~ time + Age + sex, data = CD4)
##
## Residuals:
## Min 1Q Median 3Q Max
## -531.2 -149.0 -18.7 123.7 748.3
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 479.846 71.738 6.69 0.00000000013 ***
## time 17.224 2.703 6.37 0.00000000083 ***
## Age 0.857 1.529 0.56 0.58
## sex -13.929 55.412 -0.25 0.80
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 217 on 263 degrees of freedom
## (42 observations deleted due to missingness)
## Multiple R-squared: 0.135, Adjusted R-squared: 0.125
## F-statistic: 13.7 on 3 and 263 DF, p-value: 0.000000026
#modelo asumiendo dependencia
mod_gee <- geeglm(cd4 ~ time + Age + sex , # y~ covariables
data = CD4, # datos
id = CD4$id, # id= identificador del grupo al que corresponden las medidas repetidas
family = gaussian, # modelo lineal (existen otras: poisson,gaussian, etc)
corstr = "unstructured") # tipo de correlación ("exchangeable", "independence" , "ar1", etc)
tidy(mod_gee, conf.int = T)
## # A tibble: 4 × 7
## term estimate std.error statistic p.value conf.low conf.high
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 (Intercept) 481. 63.6 57.3 3.77e-14 357. 606.
## 2 time 17.6 2.81 39.2 3.78e-10 12.1 23.1
## 3 Age 0.734 1.43 0.265 6.07e- 1 -2.06 3.53
## 4 sex -13.1 48.6 0.0722 7.88e- 1 -108. 82.3
¿COMO SE INTERPRETAN LOS COEFICIENTES? ¿QUE PASA CON LOS BETA Y LOS ERRORES STANDARD EN UNO Y OTRO MODELO? VER EFECTO DEL TIEMPO: ¿CUAL ES LA MAGNITUD DE CAMBIO DE CD4 POR MES? Y EL DE GRUPOS ¿EXISTE DIFERENCIA BASAL DE CD4 POR SEXO?
El número de cd4 en tiempo cero cuando el sexo es femenino es 581.4 El número de cd4 aumenta en 17.6 por cada mes El número de cd4 no se modifica en forma estadísticamente significativa según la edad El número de cd4 basal no se modifica en forma estadísticamente significativa según el sexo Los errores disminuyen en el modelo GEE
#GRÁFICO DE PREDICCIÓN LINEAL DEL EFECTO DEL TIEMPO SOBRE CAMBIO DE CD4 CON IC95%
#Observados
ggplot(CD4, aes(time,cd4 )) + geom_point()+
geom_smooth(method = "lm", se = TRUE) +
theme_classic()+
labs(x = "Tiempo", y = "CD4", title = "Efectos marginales")
## `geom_smooth()` using formula = 'y ~ x'
## Warning: Removed 42 rows containing non-finite outside the scale range
## (`stat_smooth()`).
## Warning: Removed 42 rows containing missing values or values outside the scale range
## (`geom_point()`).
#GRÁFICO DE PREDICCIÓN LINEAL DEL EFECTO DEL TIEMPO SOBRE CAMBIO DE CD4
# Comparo predichos
datospredichosGEE <- data.frame(na.omit(CD4), ## sacamos datos faltantes
# o no podemos graficar
phat = fitted(mod_gee))
ggplot(datospredichosGEE, aes(time, phat)) +
geom_smooth(method = "lm") +
ylim (20, 1400) +
theme_classic()+
labs(x = "Tiempo", y = "CD4", title = "Efectos marginales")
## `geom_smooth()` using formula = 'y ~ x'
Evaluacion entre sexo y tiempo
mod_gee2 <- geeglm(cd4 ~ time*sex ,
data = CD4,
id = CD4$id,
corstr = "unstructured")
tidy(mod_gee2, conf.int = T)
## # A tibble: 4 × 7
## term estimate std.error statistic p.value conf.low conf.high
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 (Intercept) 486. 75.8 41.2 1.40e-10 338. 635.
## 2 time 20.7 9.95 4.32 3.78e- 2 1.17 40.2
## 3 sex 10.1 78.5 0.0166 8.98e- 1 -144. 164.
## 4 time:sex -3.24 10.4 0.0973 7.55e- 1 -23.6 17.1
¿EL CAMBIO DE CD4 EN EL TIEMPO ES DIFERENTE EN HOMBRES Y MUJERES? no
#Evaluación del efecto del tratamiento
#Analisis exploratorio: media de CD4 en cada medición por grupo de tratamiento
CD4_WIDE %>% group_by(Treatment) %>% summarise(T0 = mean(Tiempo0, na.rm = T),
T6 = mean(Tiempo6, na.rm = T),
T12 = mean(Tiempo12, na.rm = T))
## # A tibble: 2 × 4
## Treatment T0 T6 T12
## <dbl> <dbl> <dbl> <dbl>
## 1 0 535. 601. 648.
## 2 1 468. 631. 730.
mod_gee2 <- geeglm(cd4 ~ time +Treatment+ sex ,
data = CD4,
id = CD4$id,
corstr = "unstructured")
tidy(mod_gee2, conf.int = T)
## # A tibble: 4 × 7
## term estimate std.error statistic p.value conf.low conf.high
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 (Intercept) 495. 54.1 83.7 0 389. 601.
## 2 time 17.7 2.84 39.0 4.15e-10 12.2 23.3
## 3 Treatment 11.9 26.2 0.206 6.50e- 1 -39.5 63.3
## 4 sex -7.77 48.7 0.0254 8.73e- 1 -103. 87.7
¿MIRANDO LOS VALORES A LOS 0 6 Y 12 MESES, LE PARECE QUE EXISTEN DIFERENCIAS ENTRE LOS GRUPOS? Parecería que sí
MODELO MARGINAL (PROMEDIO POBLACIONAL) DE REGRESION LINEAL ASUMIENDO DEPENDENCIA
VER EFECTO DEL TIEMPO (MAGNITUD DEL CAMBIO POR MES) Y EL DEL TRATAMIENTO (¿EXISTE DIFERENCIA BASAL DE CD4 ENTRE LOS GRUPOS?)No existe diferencia basal en los grupos de tratamiento
#gráfico de trayectoria de tratados y no tratados
mod_gee4 <- geeglm(cd4 ~ time*Treatment + sex ,
data = CD4,
id = CD4$id,
corstr = "unstructured")
tidy(mod_gee4, conf.int = T)
## # A tibble: 5 × 7
## term estimate std.error statistic p.value conf.low conf.high
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 (Intercept) 544. 53.5 104. 0 440. 649.
## 2 time 9.42 3.65 6.64 0.00998 2.25 16.6
## 3 Treatment -63.6 37.9 2.81 0.0935 -138. 10.7
## 4 sex -7.06 49.1 0.0207 0.886 -103. 89.1
## 5 time:Treatment 13.2 5.29 6.23 0.0126 2.83 23.6
#GRAFICO DEL CAMBIO DE CD4 en tratados y no tratados
predichosTTO <- data.frame(na.omit(CD4), TTO = fitted(mod_gee4))
ggplot(predichosTTO, aes(time, TTO, colour = factor(Treatment))) +
geom_smooth(method = "lm") + theme_classic()+ ggtitle("Efectos marginales") +
labs(x = "Tiempo", y = "CD4")
## `geom_smooth()` using formula = 'y ~ x'
## `geom_smooth()` using formula = 'y ~ x'
¿EL CAMBIO DE CD4 EN EL TIEMPO ES DIFERENTE EN LOS GRUPOS DE TRATAMIENTO? ¿CUÁL ES LA MAGNITUD DE CAMBIO POR MES EN CADA GRUPO? ¿CUÁL ES EL VALOR ESPERADO DE CD4 AL AÑO EN CADA GRUPO?
Es diferente. En el grupo sin tratamiento, los cd4 aumentan 9.42 cd4 por mes. en el grupo con tratamiento, aumentan 13.19+9.42=22.61