Capítulo 7 Práctica 3
7.1 Preparamos todo
# setwd("C:/your/dir")
load('../../data/Framingham2.RData')
load('../../data/Datos2.RData')
summary(Framingham) # y comprobamos lo cargado## OBS SEX CHD AGE
## Min. : 1.0 Hombre:669 Sin evidencia de CHD:1095 Min. :45.00
## 1st Qu.: 352.2 Mujer :737 Casos prevalentes : 43 1st Qu.:48.00
## Median : 703.5 Casos incidentes 4 : 39 Median :52.00
## Mean : 703.6 Casos incidentes 6 : 37 Mean :52.44
## 3rd Qu.:1054.8 Casos incidentes 7 : 31 3rd Qu.:56.00
## Max. :1406.0 Casos incidentes 10 : 30 Max. :62.00
## (Other) : 131
## SBP SBP10 DBP CHOL FRW
## Min. : 90.0 Min. : 94 Min. : 50.00 Min. : 96.0 Min. : 52.0
## 1st Qu.:130.0 1st Qu.:130 1st Qu.: 80.00 1st Qu.:200.0 1st Qu.: 94.0
## Median :142.0 Median :145 Median : 90.00 Median :230.0 Median :103.0
## Mean :148.1 Mean :148 Mean : 90.15 Mean :234.6 Mean :105.5
## 3rd Qu.:160.0 3rd Qu.:160 3rd Qu.: 98.00 3rd Qu.:264.0 3rd Qu.:114.5
## Max. :300.0 Max. :264 Max. :160.00 Max. :430.0 Max. :222.0
## NA's :635 NA's :11
## CIG YRS_CHD YRS_DTH
## Min. : 0.000 Min. : 0.00 Min. : 1.00
## 1st Qu.: 0.000 1st Qu.: 9.00 1st Qu.:18.00
## Median : 0.000 Median :16.00 Median :18.00
## Mean : 8.041 Mean :13.32 Mean :16.19
## 3rd Qu.:20.000 3rd Qu.:18.00 3rd Qu.:18.00
## Max. :60.000 Max. :18.00 Max. :18.00
## NA's :1 NA's :43
## DEATH CAUSE AGE.MONTH
## Vivo en el examen 10 :1055 Vivo en el examen 10:1056 Min. :540.0
## Fallecido en el examen 10: 60 Cáncer : 79 1st Qu.:576.0
## Fallecido en el examen 9 : 54 CHD no súbita : 60 Median :624.0
## Fallecido en el examen 6 : 49 Otras CHD : 57 Mean :629.3
## Fallecido en el examen 8 : 48 Otras enfermedades : 56 3rd Qu.:672.0
## Fallecido en el examen 5 : 38 (Other) : 79 Max. :744.0
## (Other) : 102 NA's : 19
## Genero AGE.CAT AGE.CAT2 CHD2 CIG2
## Femenino :737 50 a 54 :445 (50, 56]:523 Caso : 311 Fuma :633
## Masculino:669 54 a 59 :398 (56, 62]:343 No caso:1095 No fuma:772
## 60 o más :117 [45, 50]:540 NA's : 1
## Menos de 50:446
##
##
##
## SBP2 DBP2 CHOL2 FRW2
## Hipertenso:825 Hipertenso:717 Hipercolesterolemia:1081 Normopeso:1150
## Normotenso:581 Normotenso:689 Normocolesterolemia: 325 Sobrepeso: 245
## NA's : 11
##
##
##
## 7.2 Cálculo de porcentajes
genPercentage <- table(Framingham$SEX) *100/nrow(Framingham)
genPercentage##
## Hombre Mujer
## 47.58179 52.418217.3 Visualización de datos
7.3.1 Objetivo
Gráficas y datos numéricos para facilitar el análisis ### Framingham por sexo (cualitativa discreta) Por medio de un pie chart
pie(genPercentage, main="Distribución por sexo")
O por un gráfico de barras:
barplot(genPercentage, main="Distribución por sexo")
O por un histograma:
hist(Framingham$AGE, main="Distribución por sexo")
7.3.2 Framingham por edad (cuantitativa discreta)
boxplot(Framingham$AGE, main="Distribución de edades en el estudio") 
# tercer cuartil = percentil 757.3.3 Framingham, hábito de fumar y nivel de colesterol (cualitativa ordinal y cuantitativa continua respectivamente)
Visualizamos las tablas:
# CIG2
summary(Framingham$CIG2)*100/nrow(Framingham)## Fuma No fuma NA's
## 45.02133713 54.90753912 0.07112376# CHOL
summary(Framingham$CHOL)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 96.0 200.0 230.0 234.6 264.0 430.0La gráfica sobre el habito de fumar:
plot(Framingham$CIG2, main="Hábito de fumar en el estudio de Framingham")
La gráfica sobre el colesterol
boxplot(Framingham$CHOL, main="Nivel de colesterol en el estudio de Framingham")
7.3.4 Análisis simultáneo de varias variables
7.3.4.1 Sexo y caso de enfermedad coronaria (ambas cualitativas nominales)
chd2VsSex <- table(Framingham$CHD2, Framingham$SEX)
rbind(chd2VsSex, colSums(chd2VsSex))## Hombre Mujer
## Caso 190 121
## No caso 479 616
## 669 737# In percentage
chd2VsSexPerc <- prop.table(chd2VsSex, margin = 2)*100
rbind(chd2VsSexPerc, colSums(chd2VsSexPerc))## Hombre Mujer
## Caso 28.4006 16.41791
## No caso 71.5994 83.58209
## 100.0000 100.00000chd2VsSex <- table(Framingham$CHD2, Framingham$SEX)
barplot(chd2VsSex, beside = T) 
# var dep CHD2 y var de fila/indepe 7.3.4.2 Relación entre casos (cualitativa discreta) y la variable dependiente cualitativa discreta clasificación de peso relativo, FRW2
chd2VsFrw2 <- table(Framingham$CHD2, Framingham$FRW2)
rbind(chd2VsFrw2, colSums(chd2VsFrw2))## Normopeso Sobrepeso
## Caso 241 66
## No caso 909 179
## 1150 245# In percentage
chd2VsFrw2 <- prop.table(chd2VsFrw2, margin = 2)*100
rbind(chd2VsFrw2, colSums(chd2VsFrw2))## Normopeso Sobrepeso
## Caso 20.95652 26.93878
## No caso 79.04348 73.06122
## 100.00000 100.00000barplot(chd2VsFrw2, beside=T)
7.3.4.3 Variables agrupadas
table(Framingham$AGE.CAT, Framingham$CHD2) ##
## Caso No caso
## 50 a 54 95 350
## 54 a 59 110 288
## 60 o más 35 82
## Menos de 50 71 375# Ordenar "Menos de 50", que aparece en la última posición
# http://www.cookbook-r.com/Manipulating_data/Changing_the_order_of_levels_of_a_factor/
Framingham$AGE.CAT <- factor(Framingham$AGE.CAT, levels=c("Menos de 50", "50 a 54", "54 a 59", "60 o más"))
# Show stats
ageVsChd2 <- table(Framingham$CHD2, Framingham$AGE.CAT)
ageVsChd2##
## Menos de 50 50 a 54 54 a 59 60 o más
## Caso 71 95 110 35
## No caso 375 350 288 82summary(ageVsChd2) # return chisq.test## Number of cases in table: 1406
## Number of factors: 2
## Test for independence of all factors:
## Chisq = 21.27, df = 3, p-value = 9.254e-05ageVsChd2 <- prop.table(ageVsChd2, margin=2) * 100
rbind(ageVsChd2, colSums(ageVsChd2))## Menos de 50 50 a 54 54 a 59 60 o más
## Caso 15.91928 21.34831 27.63819 29.91453
## No caso 84.08072 78.65169 72.36181 70.08547
## 100.00000 100.00000 100.00000 100.00000Ahora ya podemos mostrar la relación entre franjas de edad y caso de enfermedad coronaria:
barplot(table(Framingham$AGE.CAT, Framingham$CHD2), beside=T) 
# "parece que el número de casos aumenta con la edad"Podemos expresar la tabla de contingencia usando el gráfico en mosaico:
pacman::p_load("graphics")
mosaicplot(ageVsChd2, shade = TRUE, las=2,
main = "Relación entre edad y enfermedad coronaria")7.3.4.4 Relación entre casos, var. cuantitativa discreta y var. cualitativa
# Porcentaje de cigarrillos al día en base al sexo
# sexVsCig <- table(Framingham$CIG, Framingham$SEX)
sexVsCig <- Framingham[,c("SEX", "CIG")]
sexVsCigPerc<-round(prop.table(table(sexVsCig), margin = 2)*100, 2)
rbind(sexVsCigPerc, colSums(sexVsCigPerc))## 0 1 5 10 15 20 25 30 33 35 40 45 50
## Hombre 33.16 24.14 40.32 49.09 60 76.19 84 84.31 100 75 96.55 50 100
## Mujer 66.84 75.86 59.68 50.91 40 23.81 16 15.69 0 25 3.45 50 0
## 100.00 100.00 100.00 100.00 100 100.00 100 100.00 100 100 100.00 100 100
## 60
## Hombre 100
## Mujer 0
## 100# Distribución de cigarrillos al día en hombres
# methd. 1
data <- split(Framingham$CIG, Framingham$SEX)
# typeof(data) check type list
summary(data$Hombre) # distribution in men## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.00 0.00 10.00 12.91 20.00 60.00summary(data$Mujer) #distribution in women## Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
## 0.000 0.000 0.000 3.617 5.000 45.000 1# meth. 2
# names(data) # rownames but for list
lapply(split(sexVsCig, sexVsCig$SEX), summary) # returns list (sapply does amtrix)## $Hombre
## SEX CIG
## Hombre:669 Min. : 0.00
## Mujer : 0 1st Qu.: 0.00
## Median :10.00
## Mean :12.91
## 3rd Qu.:20.00
## Max. :60.00
##
## $Mujer
## SEX CIG
## Hombre: 0 Min. : 0.000
## Mujer :737 1st Qu.: 0.000
## Median : 0.000
## Mean : 3.617
## 3rd Qu.: 5.000
## Max. :45.000
## NA's :1Y la gráfica:
boxplot(data)
Con la libreria ggplot2:
library(ggplot2)
plot <- ggplot(Framingham, aes(SEX, CIG, color=SEX)) + geom_boxplot() + stat_summary(fun=mean, geom="point", shape=10, size=4) # añadir media
plot + labs(SEX="Sexo")## Warning: Removed 1 rows containing non-finite values (`stat_boxplot()`).## Warning: Removed 1 rows containing non-finite values (`stat_summary()`).
7.3.4.5 Entre sexo y colesterol CHOL
# sexVsChd <- split(Framingham$CHD, Framingham$SEX) # groups=
# total por filas y porcentaje
sexVsChd <- table(Framingham$CHD, Framingham$SEX)
# https://stackoverflow.com/a/15866750/11688263
sexVsChd <- cbind(sexVsChd, rowSums(sexVsChd), rowSums(sexVsChd)*100 / sum(sexVsChd))
# y por columnas
sexVsChd <- rbind(sexVsChd, colSums(sexVsChd[,-c(3, 4)]), colSums(sexVsChd[,-c(3, 4)]*100/sum(sexVsChd[,-c(3,4)])))
sexVsChd[nrow(sexVsChd), c(3,4)]<-NaN # rm
sexVsChd[nrow(sexVsChd)-1, c(3,4)]<-NaN # rm
colnames(sexVsChd) <- c(colnames(sexVsChd)[1:2], "Total", "Perc")
sexVsChd## Hombre Mujer Total Perc
## Sin evidencia de CHD 479.00000 616.00000 1095 77.880512
## Casos prevalentes 26.00000 17.00000 43 3.058321
## Casos incidentes 2 20.00000 5.00000 25 1.778094
## Casos incidentes 3 15.00000 7.00000 22 1.564723
## Casos incidentes 4 25.00000 14.00000 39 2.773826
## Casos incidentes 5 15.00000 14.00000 29 2.062589
## Casos incidentes 6 26.00000 11.00000 37 2.631579
## Casos incidentes 7 17.00000 14.00000 31 2.204836
## Casos incidentes 8 15.00000 12.00000 27 1.920341
## Casos incidentes 9 16.00000 12.00000 28 1.991465
## Casos incidentes 10 15.00000 15.00000 30 2.133713
## 669.00000 737.00000 NaN NaN
## 47.58179 52.41821 NaN NaN# entiende que H, M y Total son float /!
# lapply(prop.table(sexVsChd), summary)Podemos expresar la tabla de contingencia usando el gráfico de barras:
barplot(table(Framingham$SEX, Framingham$CHD), beside = TRUE,
main = "Relación entre sexo y casos incidentes")
De nuevo, usando ggplot2:
ggplot(Framingham, aes(x=CHD, fill=SEX)) + geom_bar(position = "dodge") + theme(axis.text.x = element_text(angle=90)) + labs(title="Frec. de casos de enfermedad coronaria, según número y sexo", fill ="Sexo", x="Número de casos de enfermedad coronaria", y="Frecuencia absoluta")
7.3.4.6 Entre variable cuantitativa y cuantitativa
# variance
var(Framingham[, c("AGE", "CIG")], use="complete")## AGE CIG
## AGE 22.926380 -7.782218
## CIG -7.782218 134.372285#covariance
cov(Framingham[, c("AGE", "CIG")], use="complete")## AGE CIG
## AGE 22.926380 -7.782218
## CIG -7.782218 134.372285#corelation
cor(Framingham[, c("AGE", "CIG")], use="complete")## AGE CIG
## AGE 1.0000000 -0.1402106
## CIG -0.1402106 1.0000000plot(Framingham$CHOL~Framingham$AGE) # ~ left side is the dependent
pacman::p_load("ggpubr")
library(ggpubr)
ggscatter(Framingham, x = "AGE", y = "CHOL", add = "reg.line", conf.int = TRUE, main="Relación entre edad y colesterol", cor.coef = TRUE, xlab = "Edad (años)", ylab = "Nivel de colesterol")