PEC 1

Elena Alcántara Zamora

Sección 1. Importación, exportación y gestión de datos.

En primer lugar, antes de comenzar con los ejercicios, debemos descargar el archivo heart.csv.


**1.1. Importe los datos del archivo en un data frame llamado HeartStat y muestre los 6 primeros registros (head) y los 6 últimos (tail).**



``` r
HeartStat <- read.csv("Heart_disease_statlog.csv")
head(HeartStat, 6)
tail(HeartStat, 6)

1.2.Describa las características del conjunto de datos: nombre de las variables, tipos de datos, tamaño de la muestra (número de filas y columnas) y prospección de valores nulos. Proponga qué tipo de análisis exploratorio y gráficos serán adecuados según la naturaleza de las variables.

cat("names(): muestra los nombres de las variables.\n")
## names(): muestra los nombres de las variables.
names(HeartStat)
##  [1] "age"      "sex"      "cp"       "trestbps" "chol"     "fbs"     
##  [7] "restecg"  "thalach"  "exang"    "oldpeak"  "slope"    "ca"      
## [13] "thal"     "target"
cat("str(): muestra los tipos de datos.\n")
## str(): muestra los tipos de datos.
str(HeartStat)
## 'data.frame':    270 obs. of  14 variables:
##  $ age     : int  70 67 57 64 74 65 56 59 60 63 ...
##  $ sex     : int  1 0 1 1 0 1 1 1 1 0 ...
##  $ cp      : int  3 2 1 3 1 3 2 3 3 3 ...
##  $ trestbps: int  130 115 124 128 120 120 130 110 140 150 ...
##  $ chol    : int  322 564 261 263 269 177 256 239 293 407 ...
##  $ fbs     : int  0 0 0 0 0 0 1 0 0 0 ...
##  $ restecg : int  2 2 0 0 2 0 2 2 2 2 ...
##  $ thalach : int  109 160 141 105 121 140 142 142 170 154 ...
##  $ exang   : int  0 0 0 1 1 0 1 1 0 0 ...
##  $ oldpeak : num  2.4 1.6 0.3 0.2 0.2 0.4 0.6 1.2 1.2 4 ...
##  $ slope   : int  1 1 0 1 0 0 1 1 1 1 ...
##  $ ca      : int  3 0 0 1 1 0 1 1 2 3 ...
##  $ thal    : int  1 3 3 3 1 3 2 3 3 3 ...
##  $ target  : int  1 0 1 0 0 0 1 1 1 1 ...
cat("dim(): muestra el tamaño de la muestra.\n")
## dim(): muestra el tamaño de la muestra.
dim(HeartStat)
## [1] 270  14
cat("sum(is.na()): muestra los valores nulos.\n")
## sum(is.na()): muestra los valores nulos.
sum(is.na(HeartStat))
## [1] 0

El tipo de análisis explotatorio y los gráficos adecuados varían dependiendo del tipo de variable:

  • Cuantitativa continua como age, chol o oldpeak: media, mediana, desviación típica, histograma, diagrama de cajas…
  • Cualitativa binaria como sex, fbs o target: tabla de frecuencias y gráfico de barras.
  • Cualitativa con varias categorías como cp, slope o thal: tabla de frecuencias, gráfico de barras, tablas cruzadas…
  • Discreta como ca: tabla de frecuencias y gráfico de barras.

1.3. Supongamos que en este dataset detectamos valores ausentes en una variable. La práctica habitual de “eliminar todas las filas con un valor ausente” (complete-case analysis) puede introducir un sesgo. Explique qué riesgo puede comportar descartar sistemáticamente a pacientes que no tienen el expediente médico completo.

Eliminar sistemáticamente a los pacientes que tienen algún dato ausente puede introducir sesgos de selección ya que podríamos eliminar pacientes con datos diferentes de los que tienen el expediente completo.

Por ejemplo, muchos de los pacientes con expedientes incompletos son personas mayores, que suelen enfermar y acudir a urgencias con mucha más frecuencia. Si eliminamos a estos pacientes, podría ser que los expedientes completos pertenezcan a pacientes más sanos, introduciendo así sesgos en los estudios ya que no se representan bien a todos los pacientes. Así, perderíamos mucha información y podríamos pensar que hay menos personas enfermas que acuden a urgencias de las que hay en realidad.

Al eliminar los datos incompletos, podríamos estudiar solamente una parte de los pacientes, obteniendo así un resultado equivocado de la realidad.

2.1. Explique qué información contiene el conjunto de datos esoph y sus características más relevantes. Sugiera y justifique alguna acción de preparación o transformación de los datos (p. ej., creación de nuevas variables categóricas o agrupación de niveles de consumo de tabaco/alcohol).

En primer lugar, debemos cargar el dataset interno de R ‘esoph’:

data("esoph")
help("esoph")
## starting httpd help server ... done
str(esoph)
## 'data.frame':    88 obs. of  5 variables:
##  $ agegp    : Ord.factor w/ 6 levels "25-34"<"35-44"<..: 1 1 1 1 1 1 1 1 1 1 ...
##  $ alcgp    : Ord.factor w/ 4 levels "0-39g/day"<"40-79"<..: 1 1 1 1 2 2 2 2 3 3 ...
##  $ tobgp    : Ord.factor w/ 4 levels "0-9g/day"<"10-19"<..: 1 2 3 4 1 2 3 4 1 2 ...
##  $ ncases   : num  0 0 0 0 0 0 0 0 0 0 ...
##  $ ncontrols: num  40 10 6 5 27 7 4 7 2 1 ...
summary(esoph)
##    agegp          alcgp         tobgp        ncases         ncontrols     
##  25-34:15   0-39g/day:23   0-9g/day:24   Min.   : 0.000   Min.   : 0.000  
##  35-44:15   40-79    :23   10-19   :24   1st Qu.: 0.000   1st Qu.: 1.000  
##  45-54:16   80-119   :21   20-29   :20   Median : 1.000   Median : 4.000  
##  55-64:16   120+     :21   30+     :20   Mean   : 2.273   Mean   : 8.807  
##  65-74:15                                3rd Qu.: 4.000   3rd Qu.:10.000  
##  75+  :11                                Max.   :17.000   Max.   :60.000
dim(esoph)
## [1] 88  5
  • El comando help() abre un panel que incluye la descripción del estudio, el significado de cada variable, sus niveles y la fuente de los datos. Como vemos se trata de un estudio caso-control de cáncer de esófago en Ille-et-Vilaine, Francia. Contiene datos de 88 combinaciones de grupos de edad, alcohol y tabaco.
  • El comando str() nos muestra los tipos de datos. En este caso, vemos 5 variables y 88 observaciones.
  • El comando summary() nos muestra un resumen estadístico de las variables. En este caso, nos dará las frecuencias de cada nivel para agegp, alcgp y tobgp, además del resumen para ncases y ncontrols.
  • El comando dim() nos dará el tamaño de los datos. En este caso, 88 filas y 5 columnas.

Una posible transformación puede ser crear la proporción de casos. Es decir, creamos dos columnas nuevas: una que contenga el nº de personas de cada estrato y otra que sea la fracción de ellas con cáncer.

esoph_modificado <- esoph
esoph_modificado$ntotal <- esoph_modificado$ncases + esoph_modificado$ncontrols

esoph_modificado$prop_casos <- esoph_modificado$ncases / esoph_modificado$ntotal

Como ya hemos dicho, cada fila es un estrato, no corresponde a un paciente, y tienen tamaños muy distintos. Por tanto, con la proporción es mucho más fácil comparar estratos entre sí para que los resultados sean fiables.

2.2. Exporte el data frame esoph modificado a un archivo en formato .csv (llamado esoph_modificat.csv) especificando la ruta de destino y asegurándose de no incluir el índice de filas.

ruta <- "C:/Users/Usuario/ESCRITORIO/UNIVERSIDAD/MASTER BIOINFORMATICA/SOFTWARE PARA EL ANÁLISIS DE DATOS/RETO1/esoph_modificado.csv"
write.csv(esoph_modificado, file = ruta, row.names = FALSE)
file.exists(ruta)
## [1] TRUE

En primer lugar, definimos la ruta de destino como ‘ruta’. Tras esto y gracias al comando write.csv, exportamos el data frame modificado a la ruta de destino que hemos especificado. Con row.names = FALSE evitamos incluir el índice de las filas.

2.3. El dataset esoph recoge información sobre hábitos de consumo de alcohol y tabaco relacionados con una enfermedad grave. Cuando los algoritmos de salud utilizan datos de estilo de vida para evaluar riesgos o priorizar tratamientos, ¿cuál es el límite ético entre la prevención médica y la discriminación o estigmatización de determinados colectivos?

En mi opinión, el límite ético entre la prevención y la discriminación estaría en utilizar los datos de estilo de vida para mejorar la salud y no para culpabilizarlos o juzgarlos por sus hábitos o circunstancias personales. Debemos usar los datos para ayudar y proteger a los pacientes y no convertirlos en una herramienta de discriminación.

Por ejemplo, un algoritmo que utilice estos datos puede ayudar a establecer asociacione entre el consumo de alcohol, tabaco y el cáncer de esófago. Esto es muy importante a la hora de realizar campañas de prevención o detectar personas que podrían beneficiarse de un seguimiento cercano. Sin embargo, existiría un problema ético si se usan estos resultados para reducir la atención sanitaria o responsabilizar al paciente de su enfermedad por sus hábitos.

Sección 2. Análisis básico del conjunto de datos.

3.1. Calcule la proporción de pacientes según el tipo de dolor torácico (cp), agrupados por la presencia o ausencia de enfermedad cardíaca (target). Comente los resultados.

proporcion_cp <- prop.table(table(HeartStat$target, HeartStat$cp), 1)
proporcion_cp
##    
##              0          1          2          3
##   0 0.10000000 0.23333333 0.41333333 0.25333333
##   1 0.04166667 0.05833333 0.14166667 0.75833333
round(proporcion_cp, 3)
##    
##         0     1     2     3
##   0 0.100 0.233 0.413 0.253
##   1 0.042 0.058 0.142 0.758
rowSums(proporcion_cp)
## 0 1 
## 1 1

Para calcular la proporción que pide el ejercicio, uso el comando prop.table para obtener la tabla de frecuencias relativas a partir de las frecuencias absolutas (table) de las variables target (filas) y cp (columnas). Con el comando round() redondeamos el resultado a 3 decimales y con rowSums() se suman todas las proporciones para ver si el resultado final es igual a 1. En cuanto a los resultados que se obtienen, vemos que en los pacientes no enfermos, predomina el dolor tipo 2, mientras que en los pacientes enfermos, predomina el dolor tipo 3 de forma muy clara. Además, vemos que en los pacientes no enfermos los tipos de dolor se reparte de forma más o menos equilibrada, mientras que en los pacientes enfermos predomina el tipo 3 con un 75.8% de los casos.

3.2. Muestre la relación entre el sexo de los pacientes (sex: 1 = hombre, 0 = mujer) y la presencia de enfermedad cardíaca (target). Represente un gráfico de barras que compare estas proporciones.

En primer lugar, para ver la relación entre el sexo y la presencia de enfermedad cardiaca es necesario ver las proporciones dentro de cada sexo. Para ello, hacemos las proporciones a partir de la tabla de frecuencias absolutas, en la que las filas serán la variable sexo y las columnas, la variable target:

proporcion_sexo <- prop.table(table(HeartStat$sex, HeartStat$target), 1)

round (proporcion_sexo, 3)
##    
##         0     1
##   0 0.770 0.230
##   1 0.454 0.546
rowSums((proporcion_sexo))
## 0 1 
## 1 1

Al igual que hicimos anteriormente, redondeamos a 3 decimales y vemos si la suma de cada fila nos da como resultado 1.

Así, vemos que el 23% de las mujeres presentan enfermedad cardiaca frente a un 54,6% de hombres. Así, vemos que la proporción en hombres es mucho mayor, lo que podría sugerir una asociación entre el sexo y la presencia de enfermedad.

Para poder ver estas diferencias de forma visual, realizamos un gráfico de barras:

barplot(t(proporcion_sexo), beside = TRUE, 
        names.arg = c("Mujer", "Hombre"), 
        col = c("blue", "red"), 
        main = "Enfermedad cardíaca según sexo", 
        xlab = "Sexo", ylab="Proporción enfermedad")

3.3. Evalúe la presencia de enfermedad cardíaca según los niveles de azúcar en sangre en ayunas (fbs), contando el número absoluto de pacientes por cada categoría.

Para ver la presencia de enfermedad en función de los niveles de azúcar, lo hacemos gracias a una tabla de frecuencias absolutas:

tabla_fbs <- table(HeartStat$fbs, HeartStat$target)
tabla_fbs
##    
##       0   1
##   0 127 103
##   1  23  17

Así, podemos ver que hay muchos más pacientes que presentan una glucemia normal que con glucemia elevada. Para poder evaluar la presencia de enfermedad en función de la glucemia, deberíamos comparar la proporción de enfermos dentro de cada grupo, tal y como hemos hecho en apartados anteriores.

round(prop.table(tabla_fbs, 1), 3)
##    
##         0     1
##   0 0.552 0.448
##   1 0.575 0.425

Así, podemos ver que, a priori, no hay mucha diferencia entre tener glucemia alta y enfermedad. Prácticamente todas las proporciones son bastante similares.

3.4. Muestre en una tabla el número de pacientes con diagnóstico positivo de enfermedad cardíaca, agrupados por tipo de angina de esfuerzo (exang), que sean hombres y cuya presión arterial en reposo sea superior a la media de la muestra.

En primer lugar, debemos conocer cuál es la media de la presión arterial:

media_pa <- mean(HeartStat$trestbps)
media_pa
## [1] 131.3444

Ahora, debemos filtrar por diagnóstico positivo para enfermedad cardiaca, que sean hombres y que la presión anterial sea mayor a la media:

filtrado <- subset(HeartStat, target ==1 & sex==1 & trestbps > media_pa)
filtrado

Una vez filtrados, debemos agrupar por tipo de angina:

tabla_exang <- table(filtrado$exang)
tabla_exang
## 
##  0  1 
## 18 24

Así, vemos que de los pacientes que hemos filtrado, 18 no tienen angina de esfuerzo y 24 sí que la presentan.

3.5. En muchos modelos aplicados a la cardiología, la variable “sexo” o “etnia” mejora la capacidad predictiva del modelo, pero puede generar sesgos de atención si el modelo se ha entrenado mayoritariamente con un solo perfil de paciente (por ejemplo, hombres de mediana edad). ¿Qué precauciones deberían tomarse antes de implementar un algoritmo de ayuda al diagnóstico cardíaco en un hospital real?

Antes de implementar algoritmos para ayudar al diagnóstico, tendríamos que estar seguros de tener representación equitativa de toda la población. Es decir, si usamos la variable sexo, deberíamos incluir en el estudio tanto a hombres como mujeres; de igual forma con la variable etnia. Un modelo entrenado mayoritariamente con hombres de mediana edad puede ser preciso en promedio, pero fallar en mujeres o en personas mayores. Esto provocaría fallos en el diagnóstico de estas poblaciones no incluidas en el estudio. Además, es importante evaluar el rendimiento del modelo en cada grupo de pacientes, analizando posibles falsos negativos y falsos positivos. Por último, también es importante que se pruebe en condiciones reales antes de implantarlo para así poder revisar posibles fallos.

4.1. Obtenga las principales estadísticas descriptivas (media, desviación estándar, mediana, cuartiles y rango intercuartílico) de las variables cuantitativas continuas del dataset (age, trestbps, chol, thalach). Comente los resultados destacando si observa asimetría o valores atípicos.

Para obtener los principales estadísticos descriptivos, utilizamos lo siguiente:

variables <- c("age", "trestbps", "chol", "thalach")

sapply(HeartStat[, variables], mean)
##       age  trestbps      chol   thalach 
##  54.43333 131.34444 249.65926 149.67778
sapply(HeartStat[, variables], sd)
##       age  trestbps      chol   thalach 
##  9.109067 17.861608 51.686237 23.165717
sapply(HeartStat[, variables], median)
##      age trestbps     chol  thalach 
##     55.0    130.0    245.0    153.5
sapply(HeartStat[, variables], quantile, probs = c(0.25, 0.75))
##     age trestbps chol thalach
## 25%  48      120  213     133
## 75%  61      140  280     166
sapply(HeartStat[, variables], IQR)
##      age trestbps     chol  thalach 
##       13       20       67       33

Así, podemos obtener todos los parámetros estadísticos que se piden de las variables cuantitativas.

Para detectar asimetría y valores atípicos, nos fijamos en la media y mediana, además de los cuartiles. Por ejemplo, la variable age tiene una media y mediana casi iguales, y los cuartiles Q1 y Q3 están a una distancia parecida de la mediana, por lo que esta variable presenta una distribución prácticamente simétrica. Al igual ocurre con trestbps. Sin embargo, en chol vemos una media ligeramente mayor a la mediana y un Q3 más lejos de la mediana que el Q1, por lo que podríamos decir que hay una ligera asimetría positiva. Igualmente ocurre con thalach, pero presentando esta una asimetría negativa.

Para detectar valores atípicos, sería conveniente realizar un diagrama de cajas:

boxplot(HeartStat[, variables], 
        main = "Variables cuantitativas", 
        col = c("red", "blue", "orange", "yellow"))

Así, vemos que hay valores atípicos en todas las variables menos en age, siendo estos más frecuentes en la variable chol.

4.2. Realice un diagrama de cajas que relacione la frecuencia cardíaca máxima alcanzada con la presencia o ausencia de enfermedad cardíaca. ¿Qué conclusión puede extraer?

boxplot (thalach ~ target, data = HeartStat, 
         names = c("Sin enfermedad", "Enfermo"), 
         col = c("blue", "red"), 
         main = "Frec. cardíaca máx según enfermedad", 
         xlab = "Enfermedad cardíaca", ylab = "Frec. máxima")

Gracias al boxplot, podemos ver que los valores con y sin enfermedad son bastante similares, siendo un poco mayores sin enfermedad. Así, podríamos decir que los pacientes con enfermedad cardíaca tienen valores de frecuencia cardíaca máxima menores que los pacientes sin enfermedad, aunque la diferencia no es del todo clara.

4.3. Construya una tabla de frecuencias cruzadas que relacione el resultado del electrocardiograma en reposo con la presencia de enfermedad cardíaca. A continuación, represente esta información mediante un mapa de calor (heatmap) para identificar visualmente qué patrones son los más frecuentes.

tabla_cruzada <- table(HeartStat$restecg, HeartStat$target)
tabla_cruzada
##    
##      0  1
##   0 85 46
##   1  1  1
##   2 64 73

En esta tabla de contingencia, vemos la variable restecg (resultado del electrocardiograma) en filas, mientras que la variable target se encuentra en las columnas. Para poder identificar patrones frecuentes, utilizamos los gráficos de calor.

heatmap (as.matrix(tabla_cruzada),
         col = heat.colors(10), scale = "none",
         main = "ECG según enfermedad", 
         xlab = "Enfermedad cardíaca", ylab = "ECG")

Vemos que hay muy pocos pacientes con ECG 1, por lo que es poco frecuente. Sin embargo, en ECG 0 hay más pacientes sin enfermedad cardíaca que con enfermedad. Por último, en ECG 2 vemos que hay muchos pacientes y que se reparte de forma equilibrada entre enfermos y no enfermos.

4.4. Represente un histograma de la variable nivel de colesterol (chol) para los pacientes con enfermedad cardíaca. Superponga la curva de densidad empírica y la curva de densidad normal teórica. Evalúe si la variable sigue una distribución normal y comente el resultado.

En primer lugar seleccionamos los pacientes con enfermedad cardíaca.

enfermos <- subset (HeartStat, target > 0)

colesterol <- enfermos$chol

hist(colesterol, probability = TRUE, 
     main ="Colesterol en enfermos", 
     xlab = "Colesterol",
     ylab = "Frecuencia",
     col = "blue")

#Esto dibujaría la curva de densidad empírica
lines(density(colesterol), 
      col = "lightblue")

#Esta variable contiene los valores mínimos y máximos de la variable colesterol. 
x <- seq (min(colesterol), max(colesterol))

#Así genemos la curva de densidad de distribución normal teórica.
lines (x, dnorm(x, mean = mean (colesterol), sd = sd(colesterol)), 
       col = "red")

Como vemos, la curva empírica y teórica coinciden prácticamente en su totalidad, por lo que podríamos decir que la variable sigue una distribución normal, con una mayor concentración de observaciones entre unos 200-300 de colesterol.

4.5. Calcule la matriz de correlación lineal entre las variables cuantitativas del dataset y represéntela mediante un gráfico de correlación. Comente qué dos variables muestran la correlación más fuerte (positiva o negativa) e interprete su significado clínico.

Las variables cuantitativas son age, trestbps, chol, thalach y oldpeak. Creamos una variable que contenga estas 5:

variables_cuanti <- c("age", "trestbps", "chol", "thalach", "oldpeak")

#Para hacer la matriz de correlación:
matriz_cor <- cor(HeartStat[, variables_cuanti])
round (matriz_cor, 2)
##            age trestbps  chol thalach oldpeak
## age       1.00     0.27  0.22   -0.40    0.19
## trestbps  0.27     1.00  0.17   -0.04    0.22
## chol      0.22     0.17  1.00   -0.02    0.03
## thalach  -0.40    -0.04 -0.02    1.00   -0.35
## oldpeak   0.19     0.22  0.03   -0.35    1.00

El comando cor() calcula el coeficiente de Pearson, que mide la relación lineal. Aquellos valores cercanos a 1 indican una relación positiva fuerte, mientras que valores cercanos a -1 indican una relación negativa fuerte y cercanos a 0 sin relación lineal.

Para hacer el gráfico de correlación:

library(corrplot)
## corrplot 0.95 loaded
corrplot (matriz_cor, method = "color",
          title = "Matriz de correlación")

Hemos utilizado el comando corrplot de su librería, la cual hemos tenido que instalar. En cuanto a las variables con una correlación más fuerte, vemos que age y thalach son las que presentan unos círculos más grandes y oscuros, en color salmón, por lo que podemos decir que tienen correlación negativa. Esto podría indicar que, a mayor edad, menor es la frecuencia cardíaca máxima alcanzada. Esto es esperable ya que, con la edad, la frecuencia cardíaca máxima disminuye debico a una menor respuesta del corazón.