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(readxl)
library("ggplot2")
library(psych)
## 
## Adjuntando el paquete: 'psych'
## 
## The following objects are masked from 'package:ggplot2':
## 
##     %+%, alpha
library(tableone)
library(table1)
## 
## Adjuntando el paquete: 'table1'
## 
## The following objects are masked from 'package:base':
## 
##     units, units<-
library(dplyr)
options(scipen = 999, digits = 3, encoding = 'UTF-8')
library(descr)
library(expss)
## Cargando paquete requerido: maditr
## 
## To aggregate data: take(mtcars, mean_mpg = mean(mpg), by = am)
## 
## 
## Adjuntando el paquete: 'maditr'
## 
## The following objects are masked from 'package:dplyr':
## 
##     between, coalesce, first, last
## 
## The following object is masked from 'package:purrr':
## 
##     transpose
## 
## The following object is masked from 'package:readr':
## 
##     cols
## 
## 
## Use 'expss_output_rnotebook()' to display tables inside R Notebooks.
##  To return to the console output, use 'expss_output_default()'.
## 
## 
## Adjuntando el paquete: 'expss'
## 
## The following objects are masked from 'package:stringr':
## 
##     fixed, regex
## 
## The following objects are masked from 'package:dplyr':
## 
##     compute, contains, na_if, recode, vars, where
## 
## The following objects are masked from 'package:purrr':
## 
##     keep, modify, modify_if, when
## 
## The following objects are masked from 'package:tidyr':
## 
##     contains, nest
## 
## The following object is masked from 'package:ggplot2':
## 
##     vars
library(DescTools)
## 
## Adjuntando el paquete: 'DescTools'
## 
## The following object is masked from 'package:maditr':
## 
##     %like%
## 
## The following objects are masked from 'package:psych':
## 
##     AUC, ICC, SD
library(emmeans)
## Welcome to emmeans.
## Caution: You lose important information if you filter this package's results.
## See '? untidy'
library(summarytools)
## 
## Adjuntando el paquete: 'summarytools'
## 
## The following objects are masked from 'package:descr':
## 
##     descr, freq
## 
## The following objects are masked from 'package:table1':
## 
##     label, label<-
## 
## The following object is masked from 'package:tibble':
## 
##     view
library(data.table)
## 
## Adjuntando el paquete: 'data.table'
## 
## The following object is masked from 'package:DescTools':
## 
##     %like%
## 
## The following objects are masked from 'package:expss':
## 
##     copy, like
## 
## The following objects are masked from 'package:maditr':
## 
##     copy, dcast, let, melt
## 
## 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(epiR)
## Cargando paquete requerido: survival
## Package epiR 2.0.77 is loaded
## Type help(epi.about) for summary information
## Type browseVignettes(package = 'epiR') to learn how to use epiR for applied epidemiological analyses
library(sjPlot)
library(survival)
library(survminer)
## Cargando paquete requerido: ggpubr
## 
## Adjuntando el paquete: 'ggpubr'
## 
## The following object is masked from 'package:expss':
## 
##     compare_means
## 
## 
## Adjuntando el paquete: 'survminer'
## 
## The following object is masked from 'package:survival':
## 
##     myeloma
library(pwr)
library(ggpubr)
library(lmtest)
## Warning: package 'lmtest' was built under R version 4.4.3
## Cargando paquete requerido: zoo
## 
## Adjuntando el paquete: 'zoo'
## 
## The following objects are masked from 'package:data.table':
## 
##     yearmon, yearqtr
## 
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
library("openxlsx")
## Warning: package 'openxlsx' was built under R version 4.4.3
library("car")
## Cargando paquete requerido: carData
## 
## Adjuntando el paquete: 'car'
## 
## The following object is masked from 'package:DescTools':
## 
##     Recode
## 
## The following object is masked from 'package:expss':
## 
##     recode
## 
## The following object is masked from 'package:psych':
## 
##     logit
## 
## The following object is masked from 'package:dplyr':
## 
##     recode
## 
## The following object is masked from 'package:purrr':
## 
##     some
library("lmtest")
library("MASS")
## 
## Adjuntando el paquete: 'MASS'
## 
## The following object is masked from 'package:dplyr':
## 
##     select
vida <- read_excel("C:/Users/Administrador/Downloads/vida.xlsx", sheet = 1)
#nueva variable vida$genero1 tipo factor femenino=1 masc=0 
vida$genero1 <- factor(vida$genero, levels = c(0,1), labels = c("masculino", "femenino"))
#nuevo dataset solo de mujeres con nueva variable Comorbilidades tipo factor
bida=vida
bida <- bida %>% filter (genero1 %in% c("femenino"))
bida$Comorbilidades <- factor(bida$comorb, levels = c(0,1), labels = c("No", "Si"))
#tabla 1 con los datos de la muestra. descriptivo
table1(~ edad + Comorbilidades + sueño + calidad,data=bida)
Overall
(N=87)
edad
Mean (SD) 53.9 (10.5)
Median [Min, Max] 55.0 [35.0, 75.0]
Comorbilidades
No 41 (47.1%)
Si 46 (52.9%)
sueño
Mean (SD) 20.9 (2.79)
Median [Min, Max] 21.0 [15.0, 29.0]
calidad
Mean (SD) 63.8 (9.49)
Median [Min, Max] 65.0 [40.0, 85.0] 
#Evaluación de normalidad de las variables edad y calidad. Se acepta una distribución normal de ambas

describe(bida$sueño)
##    vars  n mean   sd median trimmed  mad min max range skew kurtosis  se
## X1    1 87 20.9 2.79     21    20.9 2.97  15  29    14 0.16    -0.22 0.3
hist(bida$sueño)

qqnorm(bida$sueño)
qqline(bida$sueño)

shapiro.test(bida$sueño)
## 
##  Shapiro-Wilk normality test
## 
## data:  bida$sueño
## W = 1, p-value = 0.3
describe(bida$calidad)
##    vars  n mean   sd median trimmed mad min max range  skew kurtosis   se
## X1    1 87 63.8 9.49     65      64 8.9  40  85    45 -0.26    -0.17 1.02
hist(bida$calidad)

qqnorm(bida$calidad)
qqline(bida$calidad)

shapiro.test(bida$calidad)
## 
##  Shapiro-Wilk normality test
## 
## data:  bida$calidad
## W = 1, p-value = 0.5
#Evaluación gráfica de calidad de vida y calidad de sueño. Se Realiza además coeficiente de correlación de Pearson
ggplot(bida, aes(x = sueño, y = calidad)) +
  geom_point() +
  geom_smooth(method = "lm", se = FALSE, color = "lightblue") +
  labs(title = "Relación entre Calidad de sueño y calidad de vida", x = "Calidad de sueño", y = "Calidad de vida") +
  theme_minimal()
## `geom_smooth()` using formula = 'y ~ x'

cor.test(bida$calidad, bida$sueño)
## 
##  Pearson's product-moment correlation
## 
## data:  bida$calidad and bida$sueño
## t = 9, df = 85, p-value = 0.0000000000001
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.560 0.786
## sample estimates:
##  cor 
## 0.69
#Modelo regresion lineal simple calidad de vida calidad de sueño

modelo1 <- lm(calidad~sueño, data = bida)
confint(modelo1)
##             2.5 % 97.5 %
## (Intercept)  3.45  25.90
## sueño        1.82   2.88
tab_model(modelo1)
  calidad
Predictors Estimates CI p
(Intercept) 14.68 3.45 – 25.90 0.011
sueño 2.35 1.82 – 2.88 <0.001
Observations 87
R2 / R2 adjusted 0.476 / 0.469 
#generacion de valores predichos
predichos <- fitted.values(modelo1)

#generación de residuos crudos y estandarizados
res <- residuals(modelo1)
#generación de residuos estandarizados:
standres <-rstandard(modelo1)
#gráfico de linealidad


plot(y=bida$calidad, x= bida$sueño)
abline(modelo1)

plot(y=modelo1$residuals,x=modelo1$fitted.values)
abline(h=0)

#evaluacion de normalidad

densityPlot(modelo1$residuals)

hist(modelo1$residuals)

boxplot(modelo1$residuals)

qqPlot(standres) #se utilizan residuos estandarizados

## [1]  6 79
shapiro.test(modelo1$residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  modelo1$residuals
## W = 1, p-value = 0.3
#homocedasticidad


plot(y=modelo1$residuals, x=modelo1$fitted.values)
abline(h=0)

bptest(modelo1$residuals ~ modelo1$fitted.values, data=bida)
## 
##  studentized Breusch-Pagan test
## 
## data:  modelo1$residuals ~ modelo1$fitted.values
## BP = 2, df = 1, p-value = 0.1
#identificación de outliers

standres <-rstandard(modelo1)
predichos <- modelo1$fitted.values
plot(y=standres, x=predichos)
abline(h=c(-3, 3))

standres[standres>3]
## named numeric(0)
#identificacion de valores influyentes
leverage <- as.data.frame(hatvalues(modelo1))
cooksd <- cooks.distance(modelo1)
cooksd[cooksd>(4/87)]#n es el tamaño muestral
##      1      4      6     48     56     60     79 
## 0.0679 0.0679 0.0643 0.0521 0.0517 0.0910 0.0490
#dubin watson independeencia

dwtest(modelo1)
## 
##  Durbin-Watson test
## 
## data:  modelo1
## DW = 2, p-value = 0.04
## alternative hypothesis: true autocorrelation is greater than 0
#modelo2. Puntaje de calidad de vida y género
modelo2 <- lm(calidad~genero1, data = vida)
confint(modelo2)
##                 2.5 % 97.5 %
## (Intercept)     63.04  67.42
## genero1femenino -4.33   1.53
tab_model(modelo2)
  calidad
Predictors Estimates CI p
(Intercept) 65.23 63.04 – 67.42 <0.001
genero1 [femenino] -1.40 -4.33 – 1.53 0.345
Observations 156
R2 / R2 adjusted 0.006 / -0.001 
#modelo3. Puntaje de calidad de vida y edad
modelo3 <- lm(calidad~edad, data = bida)
confint(modelo3)
##              2.5 %  97.5 %
## (Intercept) 66.916 87.4509
## edad        -0.435 -0.0608
tab_model(modelo3)
  calidad
Predictors Estimates CI p
(Intercept) 77.18 66.92 – 87.45 <0.001
edad -0.25 -0.43 – -0.06 0.010
Observations 87
R2 / R2 adjusted 0.075 / 0.065 
#modelo4.Puntaje de calidad de vida y presencia de comorbilidades
modelo4 <- lm(calidad~Comorbilidades, data = bida)
confint(modelo4)
##                  2.5 % 97.5 %
## (Intercept)       64.5  70.10
## ComorbilidadesSi -10.4  -2.78
tab_model(modelo4)
  calidad
Predictors Estimates CI p
(Intercept) 67.32 64.54 – 70.10 <0.001
Comorbilidades [Si] -6.60 -10.42 – -2.78 0.001
Observations 87
R2 / R2 adjusted 0.122 / 0.111