EL CASO BINARIO

Regresión logística (curva ROC)

Departamento de Matemáticas y Estadística, Universidad del Norte (Barranquilla, Colombia)

hllinas@uninorte.edu.co

Biographical sketch

14/08/21

Abstract

La teoría mencionada puede revisarse en el volumen 8 de mis notas de clase que aparecen en el siguiente documento: 2.2. Regresión logística y en la referencia: LLinás (2006). En Rpubs:: toc se pueden ver otros documentos de posible interés.

1 Librerías

library(aplore3)     #Base de datos para los ejemplos
library(lsm)         #Base de datos para ejemplos y estimaciones del Log-verosimilitud
library(tidyverse)   #Incluye a dplyr y ggplot2

2 Introducción

Los métodos de regresión se han convertido en un componente integral de cualquier análisis de datos preocupado por describir la relación entre una variable de respuesta y una o variables más explicativas. Muy a menudo, la variable de resultado es discreta, tomando un valor de dos o más valores posibles. El modelo de regresión logística es el más modelo de regresión de mayor uso frecuente para el análisis de estos datos.

En el documento Rpbus :: Modelos lineales generalizados se explicó que estos modelos hacen parte de los modelos lineales generalizados. Allí se pueden ver los detalles correspondientes. Para conocer con profundidad estos modelos, es importante estudiar los siguientes cuatro tipos de modelos:

• Modelo de Bernoulli.

• Modelo completo.

• Modelo nulo.

• Modelo saturado.

En el documento Rpbus :: Modelos completo, nulo y saturado se describieron sus propiedades, con los ejemplos correspondientes y en los documentos Rpbus :: Modelos logísticos (estimaciones), Rpbus :: Modelos logísticos (intervalos de confianza) y Rpbus :: Modelos logísticos (pruebas de hipótesis) se detalló todo lo relacionado con las estimaciones e intervalos de confianza para los parámetros logísticos, así como diferentes pruebas de comparación entre modelos. En este documento, se utilizarán las notaciones utilizados allá, así como los resultados encontrados en los ejemplos aplicados en esos documentos. A pesar de ello, se hará una breve descripción de los cuatro modelos mencionados arriba, los cuales serán base para la teoría que se explicará posteriormente.

3 Datasets

Para las aplicaciones, se utilizaron bases de datos de las librerías aplore3 (creado por Braglia, 2016) y lsm (creado por LLinás, Fábregas y Villalba, 2020).

library(aplore3) 
library(lsm)

Ambos paquetes incluyen, de manera “no oficial”, todos los conjuntos de datos utilizados en el texto de Hosmer, Lemeshow and Sturdivant (2013). En este link o en este otro se encuentran los nombres de los datasets con los respectivos detalles:

• Descripción (description).

• Uso (usage).

• Formato (format).

• Fuente (source) o Referencias (References).

Se resalta el hecho que lsm es un paquete que contiene:

1. Otros datasets que son de autoría propia.

2. La función lsm(), la cual nos permite estimar, entre otros, el logaritmo de la función de verosimilitud de los modelos completos, nulo, saturado y logístico. La teoría relacionada con los tres primeros modelos se van a explicar más adelante dentro de este documento. La que se relaciona con lo smodelos logísticos, en otros documentos.

Para los ejemplos, se utilizará chdage, el cual contiene datos recogidos con el fin de estudiar si la edad es un factor influyente en la presencia o no de enfermedades coronarias (CHD). Es un data frame con 100 observaciones. En aplore3 hay cuatro variables y en lsm, tres. A continuación, se describen cada una (entre paréntesis, el nombre de la variable en lsm):

1. id (ID): código de identificación (1-100).

2. chd (CHD): Presencia of CHD. En aplore3, los niveles son No, Si. En lsm son 0 (=No) y 1 (=Si).

3. age (AGE): edad (en años) de los participantes.

4. agegrp: edad agrupada (1: 20-39, 2: 30-34, 3: 35-39, 4: 40-44, 5: 45-49, 6: 50-54, 7: 55-59, 8: 60-69).

Cargamos el data frame chdage de la librería aplore3:

chdage <- aplore3::chdage
attach(chdage)

Las primeras 10 observaciones son:

ID	age	agegrp	chd
1	20	20-39	No
2	23	20-39	No
3	24	20-39	No
4	25	20-39	No
5	25	20-39	Yes
6	26	20-39	No
7	26	20-39	No
8	28	20-39	No
9	28	20-39	No
10	29	20-39	No

Observe que, en el datasets aplore3::chdage, la variable chd es binaria, pero no numérica. Al convertirla en integer, se codifica como 1=No y 2=Si. Por esta razón, esta variable se puede codificar como 0=No y 1=Si, de la siguiente manera:

CHD <- as.integer(chd)-1

4 Modelo de Bernoulli

La variable de interés \(Y\) es de Bernoulli. En símbolo, \(Y \sim {\cal B}(1,p)\), siendo

\[p:=E(Y)=P(Y=1)\]

la probabilidad de que ocurra \(Y\).

Haciendo \(n\) observaciones independientes de \(Y\), se obtiene la muestra \(Y=(Y_1, \ldots, Y_n)\) con los datos \(y_i\in\{0,1\}\), para todo \(i=1, \cdots ,n\), donde \(y_i\) es un posible valor de \(Y_i\), las cuales son independientes entre sí.

Se llega a un modelo estadístico de Bernoulli:

\[Y_i= p_i + e_i \quad\sim\quad {\cal B} (1,p_i), \quad i=1,\cdots,n\]

Fijando \(y=(y_1, \cdots ,y_n)^T\) obtenemos la función de verosimilitud en el parámetro \(p=(p_1, \cdots, p_n)^T\):

\[L(p)= \prod^n_{i=1}[p_i^{y_i}(1-p_i)^{1-y_i}]\]

Theorem 4.1 (Log-verosimilitud)

En un modelo de Bernoulli, el logaritmo de la función de máxima verosimilitud será:

\[\begin{equation} {\cal L}(p):= \ln L(p)= \sum^n_{i=1}[y_i\ln p_i + (1-y_i)\ln (1-p_i)] \tag{4.1} \end{equation}\]

Remark.

Como \(0 \le f(y,p) \le 1\), de la expresión (4.1), se tiene que

\[-\infty \le {\cal L}(p) \le 0\]

Hay varias situaciones que se pueden presentar en un modelo de Bernoulli. Se dice que éste se puede identificar como alguno de los siguientes modelos: completo, nulo o saturado.

5 Modelo completo

Definition 5.1

El modelo completo es caracterizado por el supuesto de que todos \(p_i\), \(i=1, \cdots ,n\) son considerados como parámetros.

Remark.

El siguiente teorema describe las estimaciones en este modelo:

Theorem 5.1 (Estimaciones en el modelo completo)

En el modelo completo, las ML-estimaciones de \(p_i\) son \(\hat{p}_i=Y_i\) con valores \(\hat{p}_i=y_i\),\(i=1, \cdots ,n\). Además, la estimación de la función de verosimilitud, de su logaritmo y de la llamada desviación tienen los siguientes valores:

\[{\cal L}(y) = 0, \quad L(y)= 1 \quad \mbox{y}\quad -2\,{\cal L}(y)=0\]

Example 5.1

Para los datos del archivo chdage, en el modelo completo, se tiene que \({\cal L}(y) = -2\; {\cal L}(y)= 0\) y \(L(y)=1\). En R, se pueden verificar así:

Y <- CHD
a <- ifelse(Y==0,log(1-Y),log(Y))
LogCompleto <- sum(a)
L_Completo <- exp(LogCompleto)
DevLogCompleto <- -2*LogCompleto

El paquete lsm calcula directamente el valor de \({\cal L}(y)\):

datos <- lsm::chdage
attach(datos)
modelo <- lsm(CHD~AGE, data=datos)
DevLogCompleto <- modelo$Log_Lik_Complete
DevLogCompleto

## [1] 0

6 Modelo nulo

Definition 6.1

El modelo nulo es caracterizado por el supuesto de que todos los \(p_i\), \(i=1, \cdots ,n\) son considerados iguales; es decir, se tiene un solo parámetro \(p=p_i, i=1, \cdots ,n\).

Theorem 6.1 (Log-verosimilitud en el modelo nulo)

En este caso, (4.1) será:

\[\begin{eqnarray} {\cal L}(p)&=& n[\overline{y}\ln p + (1- \overline{y})\ln (1-p)] \tag{6.1} \end{eqnarray}\]

El siguiente teorema describe las estimaciones en este modelo:

Theorem 6.2 (Estimaciones en el modelo nulo)

En el modelo nulo, la ML-estimación de \(p\) es \(\hat{p}=\overline{Y}\) con valor \(\hat{p}=\overline{y}\). Además, \[{\cal L}(\overline{y})<0 \qquad \mbox{si y sólo si}\qquad 0 < \overline{y} < 1\]

Remark.

Algunas propiedades de la gráfica de \({\cal L}(p)\) son:

1. Tiene un único mínimo en \(p=\frac{1}{2}\).

2. Teniendo en cuenta que el dominio de \({\cal L}(p)\) es el intervalo abierto \((0,1)\), la gráfica de \({\cal L}(p)\):

   1. Es decreciente en \((0,\frac{1}{2})\) y creciente en \((\frac{1}{2},1)\).

   2. No tiene puntos de inflexión en \((0,1)\).

   3. Es cóncava hacia arriba.

   4. El punto mínimo es \(\left(\frac{1}{2},-n\log 2\right)\).

El siguiente teorema indica la forma de la gráfica de \({\cal L}(p)\), para \(0<p<1\).

Theorem 6.3 (Gráfica en el modelo nulo)

La gráfica de la función \({\cal L}(p\)) en el modelo nulo es:

n <- 100
p <- seq(0.0001, 1, 0.00005)
a <- n*(p*log(p)+(1-p)*log(1-p))
D <- data.frame(p,a)

ggplot(D, aes(y = a, x = p)) +
geom_point(aes(color=a), alpha = 1.9) +
labs(x="Probabilidad p", y="Log L(p)", fill= "") + 
ylim(-n*log(2), 0)+
facet_wrap(. ~ "Gráfica del Log L(p) en el modelo nulo con n=100") +    
theme_bw(base_size = 12) +
theme(legend.position = "none")

Example 6.1

Para los datos del archivo chdage, en el modelo nulo:

1. \(\widehat{p} \;=\; \overline{y}\;=\; 43/100 \;=\; 0.43\).

2. \({\cal L}(\overline{y})\;=\; {\cal L}(0.43) \;= \; -68.33\).

3. \(-2{\cal L}(\overline{y}) \;=\; 136.66\).

4. El punto mínimo es \(\left(\frac{1}{2},-69.31\right)\).

En R, se puede verificar así:

n <- length(Y)
YBarra <- mean(Y)
LogNulo <- n*(YBarra*log(YBarra)+(1-YBarra)*log(1-YBarra))
DevNulo <- -2*LogNulo
minimo <- -n*log(2)

El paquete lsm calcula directamente el valor de \({\cal L}(\overline{y})\):

datos <- lsm::chdage
attach(datos)
modelo <- lsm(CHD~AGE, data=datos)
DevNulo <- modelo$Log_Lik_Null
DevNulo 

## [1] -68.33149

7 Modelo saturado

El modelo saturado está caracterizado por dos supuestos.

Hypothesis 7.1 (Supuesto 1 en el modelo saturado)

Se supone que:

1. Se tienen \(K\) variables explicativas \(X_1, \cdots, X_K\) (algunas pueden ser numéricas y otras categóricas) con valores \(x_{i1}, \cdots, x_{iK}\) para \(i=1, \cdots, n\) (fijadas u observadas por el estadístico, según sean variables determiní}sticas o aleatorias).

2. Entre las \(n\) kuplas \((x_{i1}, \cdots, x_{iK})\) de los valores de la variable explicativa \(X\) haya \(J\) kuplas diferentes, definiendo las \(J\) poblaciones. Por tanto, \(J \le n\).

Remark.

Para cada población \(j=1, \cdots ,J\) se denota:

• El número de observaciones \(Y_{ij}\) en cada población \(j\) por \(n_j\), siendo \(n_1+\cdots +n_J=n\);

• La suma de las \(n_j\) observaciones \(Y_{ij}\) en \(j\) por

\[Z_j:=\sum\limits_{i=1}^{n_j}Y_{ij} \quad \mbox{con valor}\quad z_j=\sum\limits_{i=1}^{n_j}y_{ij},\quad \mbox{siendo}\quad \sum\limits^J_{j=1}z_j \;= \; \sum\limits^n_{i=1}y_i\]

En la tabla 7.1 se muestra un ejemplo hipotético de un conjunto de datos que se puede organizar en \(J=4\) poblaciones.

Table 7.1: Ilustración de un cojunto de datos agrupado en \(J=4\) poblaciones

1	2	3	4	5	6	7	8	9
\(Y\)	\(X_1\)	\(X_2\)	\(X_3\)	\(X_4\)	\(X_5\)	\(j\)	\(n_j\)	\(Z_j\)
Población: Bajo, 80, Si, 170, Estrato 1
1	Bajo	80	Si	170	Estrato 1	\(j=1\)	\(n_1=3\)	\(Z_1=2\)
0	Bajo	80	Si	170	Estrato 1
1	Bajo	80	Si	170	Estrato 1
Población: Mediano, 100, Si, 150, Estrato 5
0	Mediano	100	Si	150	Estrato 5	\(j=2\)	\(n_2=4\)	\(Z_2=3\)
1	Mediano	100	Si	150	Estrato 5
1	Mediano	100	Si	150	Estrato 5
1	Mediano	100	Si	150	Estrato 5
Población: Mediano, 100, No, 180, Estrato 2
1	Mediano	100	No	180	Estrato 2	\(j=3\)	\(n_3=2\)	\(Z_3=1\)
0	Mediano	100	No	180	Estrato 2
Población: Alto, 100, No, 100, Estrato 4
0	Alto	100	No	100	Estrato 4	\(j=4\)	\(n_4=3\)	\(Z_4=2\)
1	Alto	100	No	100	Estrato 4
1	Alto	100	No	100	Estrato 4
General. \(Y\) es la variable de respuesta; \(X_1, \cdots, X_5\) son las variables explicativas; \(j\) es la població; \(n_j\) es el tamaño de la población \(j\); \(Z_j\) es el número de éxitos en la población \(j\).

Hypothesis 7.2 (Supuesto 2 en el modelo saturado)

Para mayor simplicidad en la escritura, se abreviará la j-ésima población \((x_{1j}, \cdots ,x_{Kj})\) por el símbolo \(\star\). Para cada población \(j=1, \cdots ,J\) y cada observación \(i=1,\cdots,n\) en \(j\), se supone que:

1. \((Y_{ij}|\star)\) es de Bernoulli. Es decir,

\[(Y_{ij}|\star) \sim {\cal B}(1,p_j)\]

2. Las variables \((Y_{ij}|\star)\) son independientes entre sí.

3. La esperanza y la varianza son, respectivamente,

\[p_j=P(Y_{ij}=1|\star)=E(Y_{ij}|\star), \qquad V(Y_{ij}|x_{j})=p_j(1-p_j)\]

A continuación, se oprimirá el símbolo \(\star\).

Remark.

El supuesto 2 implica:

1. Todos los \(p_{ij}\), \(i=1, \cdots ,n\) dentro de cada población \(j\) son iguales. Es decir, se tiene como parámetro el vector \(p=(p_1, \cdots ,p_J)^T.\)

2. Para cada población \(j=1, \cdots ,J\):

   • La variable \(Z_j\) es binomial. Es decir,

   \[Z_j\sim{\cal B}(n_j,p_j)\]

   • Las variables \(Z_j\) son independientes entre las poblaciones.

Theorem 7.1 (Log-verosimilitud en el modelo saturado)

En el modelo saturado, el logaritmo de la función de máxima verosimilitud será

\[\begin{eqnarray} {\cal L}(p) &=& \sum^J_{j=1}\left(\sum_{i=1}^{n_j}[{y_{ij}}\ln p_j + (1- y_{ij})\ln (1-p_j)]\right)\nonumber\\ &=& \sum^J_{j=1}[{z_j}\ln p_j + (n_j- z_j)\ln (1-p_j)] \tag{7.1} \end{eqnarray}\]

Theorem 7.2 (Estimaciones en el modelo saturado)

En el modelo saturado, las ML-estimaciones de \(p_j\) son \(\tilde{p}_j=\frac{Z_j}{n_j}\), con valores \(\tilde{p}_j=\frac{z_j}{n_j}\),\(j=1,\cdots ,J\). Además,

\[\begin{eqnarray} {\cal L}_s\;:=\;{\cal L}(\tilde p) &=& \sum^J_{j=1}n_j[\tilde{p}_j\ln \tilde{p}_j +(1-\tilde{p}_j) \ln(1-\tilde{p}_j)] \end{eqnarray}\]

También se cumple que

\[{\cal L}_s<0\quad \mbox{para}\quad 0< \tilde{p}_j <1\]

Example 7.1

Para los datos del archivo chdage, en el modelo saturado,hay \(J=43\) poblaciones y se cumple que \({\cal L}(\tilde{p})=-41.7991\), como se indica en la última fila de la Tabla 7.2:

chdage$CHD <- CHD 

chdage %>%
  group_by(age) %>%
  summarise(nj = n(),
            zj = sum(CHD)) %>%
  mutate(pj = zj/nj,
         pj = round(pj,3),
         Lp = ifelse(zj==0 | zj== nj, 0, zj*log(pj)+(nj-zj)*log(1-pj)),
         Lp = round(Lp, 4)
         ) -> saturado

Totales <- c("Total", sum(saturado$nj), sum(saturado$zj), "", sum(saturado$Lp))
Tabla <- rbind(saturado, Totales) 

Table 7.2: Estimaciones en el modelo saturado.

age	nj	zj	pj	Lp
20	1	0	0	0
23	1	0	0	0
24	1	0	0	0
25	2	1	0.5	-1.3863
:	:	:	:	:
:	:	:	:	:
63	1	1	1	0
64	2	1	0.5	-1.3863
65	1	1	1	0
69	1	1	1	0
Total	100	43		-41.7991

El paquete lsm calcula directamente los valores de \(J\) y \({\cal L}(\overline{y})\):

datos <- lsm::chdage
attach(datos)
modelo <- lsm(CHD~AGE, data=datos)
J <- modelo$Populations
devSaturado <- modelo$Log_Lik_Saturate
cbind(J, devSaturado)

##       J devSaturado
## [1,] 43   -41.79938

8 Modelo logístico

Hypothesis 8.1 (Supuesto 3: matriz de diseño)

Se hacen los supuestos 1 y 2 del modelo saturado (véase las hipótesis 7.1 y 7.2), donde adicionalmente se supone que la matriz de diseño

\[C=\left(\begin{array}{cccc} 1 & x_{11} &\cdots &x_{1K}\\ 1 & x_{21} &\cdots &x_{2K}\\ \vdots &\vdots & &\vdots\\ 1 &x_{J1} &\cdots &x_{JK}\\ \end{array}\right)\]

tiene rango completo \(Rg(C)=1+K\leq J\).

Hypothesis 8.2 (Supuesto 4: modelo logístico)

Para llegar a un modelo logístico se hace el supuesto adicional:

\[\begin{equation} \mbox{Logit}(p_j):= \ln\left(\frac{p_j}{1-p_j}\right) = \delta \;+\; \beta_1 \,x_{j1} \;+\;\cdots \;+\; \beta_K \,x_{jK} \tag{8.1} \end{equation}\]

Remark.

1. \(\alpha =(\delta,\beta_1,\ldots,\beta_K)^T\) es el vector de parámetros en el modelo.

2. Nótese que el supuesto sobre \(Rg( C)=1+K\), hace identificable al parámetro \(\alpha\).

9 Riesgo

Definition 9.1 (Riesgo)

En la práctica, la probabilidad \(p_j\) es conocida como riesgo.

Theorem 9.1 (Fórmula para el riesgo)

Sea \(g_j:=\delta \;+\; \beta_1\,x_{1j} \;+\;\cdots \;+\; \beta_K \,x_{Kj}\). Entonces, la probabilidad

\[p_j=P(Y_j=1|x_{j1}, \cdots, x_{jK})\]

de obtener un éxito en la población \(j=1, \ldots, J\), dado los valores \(x_{j1}, \cdots, x_{jK}\), viene dada por:

\[\begin{equation} p_j \;= \; \mbox{Logit}^{-1}(g_j) \;= \; \frac{e^{g_j}} {1 + e^{g_j}} \tag{9.1} \end{equation}\]

Example 9.1

En la figura 9.1 aparece la gráfica de las dos funciones siguientes: \[\begin{eqnarray*} p &=& \frac{e^{x}} {1 + e^{x}} \qquad \mbox{(Figura a, roja, con $\beta=1$)},\\ p &=&\frac{e^{-x}} {1 + e^{-x}} \qquad \mbox{(Figura b, verde, con $\beta=-1$)} \end{eqnarray*}\]

Observe que el signo de la pendiente (o sea, el valor numérico que multiplica a \(x\)) influye en el hecho si la gráfica es creciente o decreciente. La No.1 es creciente (por tener pendiente postiva) y la No.2, es decreciente (por tener pendiente negativa).

VI <- seq(-10, 10, 0.05)
Crec <- exp(VI)/(1+exp(VI))
Decr <- exp(-VI)/(1+exp(-VI))

ggplot() +
geom_point(mapping=aes(y=Crec, x = VI, color="Crec"), size=1.5) +
geom_point(mapping=aes(y=Decr, x = VI, color="Decr"),size=1.5) +  
labs(x="Variable explicativa X", y="Probabilidad de éxito p", fill= "") + 
#ylim(0, 100)+
facet_wrap(. ~ "Gráfica de dos logits con una sola variabla explicativa") +    
theme_bw(base_size = 12) +
#theme(legend.position = "none")+
  scale_color_manual(values = c("Crec" = "darkblue","Decr" = "red")) +
  labs(color = "Logits")+
scale_color_discrete(name = expression(paste("Pendiente", " ", beta, ":")), labels = c(expression(paste("(a)", " ", beta>0)), expression(paste("(b)", " ", beta<0)))) # Edit legend title and labels

Figure 9.1: Comparación de dos Logits

Theorem 9.2 (Log de la función de verosimilitud en el modelo logístico)

Sea \(g_j:=\delta \;+\; \beta_1\,x_{1j} \;+\;\cdots \;+\; \beta_K \,x_{Kj}\). Reescribiendo \({\cal L}(p)\), dada en la ecuación (7.1), el logaritmo de la función de verosimilitud se puede escribir, en función de \(\alpha\), como:

\[\begin{eqnarray} {\cal L}({\alpha}) &=& \sum^{J}_{j=1}\left[z_j\ln\left(\frac{p_j}{1-p_j}\right)+ n_j\ln(1-p_j)\right]\\ &=& \sum^{J}_{j=1}z_j\,g_j \;-\; \sum^{J}_{j=1} n_j\ln \left[1+ e^{g_j}\right] \tag{9.2} \end{eqnarray}\]

En la primera expresión se observa que el paso de \(p_j\) hacia \(\ln\left(\frac{p_j} {1-p_j}\right)\) aparece de una manera natural.

10 Odds

Definition 10.1 (Odds)

Un odds se define como la proporción entre las probabilidades de ocurrencia y no ocurrencia del evento que se relaciona con \(Y\) en la población \(j\). Es decir, es el cociente

\[O_j\;=\;\frac{p_j}{1-p_j}\]

Remark (relación entre Riesgo y Odds).

Se resalta el hecho de que los riesgos toman valores entre 0 y 1; los odds, entre 0 e infinito. . Además, observe que:

\[p_j\;=\;\frac{O_j}{1+O_j}\]

Example 10.1

Para más detalles, aclaraciones y ejemplos relacionados con este tema, puede consultarse el documento Rpbus :: Modelos logísticos (intervalos de confianza).

11 Riesgo relativo RR

Definition 11.1 (Riesgo relativo)

El riesgo relativo se define como el cociente entre el riesgo en un grupo con un factor de exposición o de riesgo (población \(i\)) y el riesgo en un grupo de referencia, que no tiene el factor de exposición (población \(j\)). Es decir, es el siguiente cociente:

\[RR(i\; \mbox{vs} \; j) = \frac{\mbox{Incidencia acumulada en la población $i$}}{\mbox{Incidencia acumulada en la población $j$}}\]

O sea, el cociente entre las siguientes probabilidades de éxitos:

\[RR(i\; \mbox{vs} \; j) \;=\; \frac{p_i}{p_j}=\frac{P(\mbox{$Y=1$ | población $i$})}{P(\mbox{$Y=1$ | población $j$})}\]

Example 11.1

Para más detalles, aclaraciones y ejemplos relacionados con este tema, puede consultarse el documento Rpbus :: Modelos logísticos (intervalos de confianza).

12 Razón odds

En estudios de cohortes el RR se estima de forma directa ya que se conoce la incidencia de la enfermedad en expuestos y en no expuestos. Por el contrario, en los estudios de casos y controles no se puede calcular la incidencia, porque la población de estudio se selecciona a partir de individuos que ya han desarrollado la enfermedad. Por esta razón, en los estudio de casos y controles se calcula la razón de odds u odds ratio (OR).

Definition 12.1 (Razón odds)

Una razón ODDS se define como el cociente entre dos odds:

\[OR(i\; \mbox{vs} \;j) \;=\; \frac{O_i}{O_j} \;=\; \frac{\frac{p_i}{1-p_i}}{\frac{p_j}{1-p_j}}\]

Theorem 12.1 (OR es exponencial de la pendiente)

Siempre se cumple que \(RR(i\; \mbox{vs} \;j)\) es un número entre 0 e infinito. Además, en un modelo de regresión logística, se cumple que:

\[\begin{eqnarray} OR(i\; \mbox{vs} \;j) \;= \; e^{\beta_1(x_{i1}-x_{j1}) \; + \; \beta_2(x_{i2}-x_{j2}) \; + \; \cdots \; + \; \beta_K(x_{iK}-x_{jK})} \end{eqnarray}\]

Cuando \(x_{ik}-x_{jk}=1\) para todo \(k=1, \ldots, K\), entonces

\[OR:=OR(i\; \mbox{vs} \;j) \;= \;e^{\beta_1\;+\;\cdots \;+\;\beta_K}\]

Es decir, no depende de \(X_1, \ldots, X_K\) y muestra el cambio proporcional en la variable de respuesta cuando las variables independientes se incrementen en \(1\) unidad.

Example 12.1

Para más detalles, aclaraciones y ejemplos relacionados con este tema, puede consultarse el documento Rpbus :: Modelos logísticos (intervalos de confianza).

13 Método de estimación

El método que se propone para calcular las ML-estimaciones en un modelo logístico es el método iterativo de Newton-Raphson. Generalmente, el método requiere:

1. Una estimación inicial para el valor que maximiza la función.

2. La función es aproximada en una vecindad de aquella estimación por un polinomio de segundo grado.

3. Entonces,la siguiente estimación se calcula como el máximo de dicho polinomio.

4. Luego, se repite el proceso, usando esta estimación como la estimación inicial.

5. De esta manera, el método genera una sucesión de estimaciones. Estas estimaciones convergen a la localización del máximo cuando la función es adecuada y/o la estimación inicial es buena.

Para más detalles, ver el teorema 8 en LLinás (2006). En R, las funciones glm() y lsm() calculan estas estimaciones.

14 Casos agrupado y no agrupado

1. Cuando se trabaja con el modelo saturado, se tiene el caso de utilizar datos agrupados.

2. Cuando se tiene el caso especial \(n_j=1\), para todo \(j\) (lo que implica que \(J=n\)) se habla de datos no agrupados.

3. La distinción entre datos agrupados y no agrupados es importante por dos razones:

   1. Algunos métodos de análisis apropiados a datos agrupados no son aplicables a datos no agrupados.

   2. Las aproximaciones asintóticas pueden estar basados en uno de estos dos casos distintos: (i) \(n\to\infty\) o (ii) \(J\to\infty\), caso que es únicamente es apropiado para datos no agrupados.

4. En la práctica:

   1. Cuando se tienen datos agrupados es importante tener en cuenta que \(J\) debe ser fijo. Por esta razón, debe tomarse como base el modelo saturado. Es decir, se empieza el análisis usando los vectores \(Z_j\), \(j=1,\cdots,J\).

   2. Si \(J\to\infty\) (por ejemplo, si \(J=n\)), entonces, en el modelo saturado no se puede considerar a \(J\) como fijo. Obsérvese que esta situación se presenta cuando se tienen datos no agrupados. En este caso, no se puede tomar como base el modelo saturado. Ahora se empezaría el análisis utilizando, de una vez, las observaciones \(Y_i\), \(i=1,\cdots, n\).

15 Ejemplo 1: Enunciado

Considere los datos del archivo chdage. Suponga que se quiere analizar un modelo de regresión logística, considerando a chd como variable dependiente y age como independiente.

1. Escriba, matemáticamente, el vector de parámetros logísticos y el de sus estimadores.

2. Escriba, matemáticamente, La probabilidad estimada de que un individuo tenga enfermedades coronarias (chd\(=1\)), cuando tiene una edad determinada (digamos, age\(=x_j\)).

3. Escriba, matemáticamente, el modelo logístico estimado.

4. Obtenga las estimaciones \(\hat{\delta}\) y \(\hat{\beta}\) de los parámetros logísticos \(\delta\) y \(\beta\), sin utilizar la función summary().

5. Obtenga las estimaciones \(\hat{\delta}\) y \(\hat{\beta}\) de los parámetros logísticos \(\delta\) y \(\beta\), utilizando la función summary().

6. Utilizando las estimaciones halladas en el inciso anterior, escriba en el modelo correspondiente .

7. Haga la gráfica del riesgo de tener enfermedades coronarias versus la edad. ¿Es directa o indirecta esta relación?

8. Haga la gráfica del logit estimado versus la edad. ¿Qué tipo de relación hay entres estas dos variables?

9. Estime el logit estimado, para un sujeto con 50 años.

10. Estime la proporción estimada de personas con presencia de enfermedades coronarias a la edad de 50.

11. Halle los errores estándares estimados \(\hat{S}_{\hat{\beta}}\) y \(\hat{S}_{\hat{\delta}}\) de \(\hat{\beta}\) y \(\hat{\delta}\), respectivamente.

12. Calcule \({\cal L}(\hat{\alpha})\), la estimación del logaritmo de la función de máxima verosimilitud en el modelo logístico.

13. Encuentre el odds estimado para individuos con edad de 50 e interprételo.

14. Haga la gráfica del odds estimado versus la edad.

15. Halle la razón odds estimada cuando el incremento la edad se incrementa en 1 año. Inteprétela.

16 Ejemplo 1: Solución

16.0.1 Solución parte (a)

Los vector de parámetros y de sus estimadores son, respectivamente,

\[\alpha =(\delta,\beta)^T, \qquad \hat{\alpha} =\left(\hat{\delta}, \hat{\beta}\right)^T\]

Aquí, \(T\) indica la transpuesta del vector.

16.0.2 Solución parte (b)

La probabilidad estimada de que un individuo tenga enfermedades coronarias (chd\(=1\)), cuando tiene una edad determinada (digamos, age\(=x_j\)), se puede escribir así:

\[\hat{p}_j = \hat{P}(\mbox{chd}=1 |\mbox{age}=x_j)\]

16.0.3 Solución parte (c)

Sabiendo que \(\hat{p}_j\) es como en el inciso anterior, el modelo estimado se puede escribir teniendo en cuenta la ecuación (8.1) o la (9.1):

\[\mbox{Logit}(\hat{p}_j):= \ln\left(\frac{\hat{p}_j}{1-\hat{p}_j}\right) = \hat{\delta} + \hat{\beta} x_j, \qquad \qquad \hat{p}_j \;= \; \frac{e^{\hat{\delta} + \hat{\beta} x_j}} {1 + e^{\hat{\delta} + \hat{\beta} x_j}} \]

16.0.4 Solución parte (d)

En R, las estimaciones de los parámetros logísticos \(\delta\) y \(\beta\) se pueden obtener con las funciones lsm() o glm():

# Con la función lsm:

datos <- lsm::chdage
attach(datos)
modelo <- lsm(CHD~AGE, data=datos)
modelo$coef     # Todos los coeficientes
modelo$coef[1]  # Intercepto
modelo$coef[2]  # Pendiente

# Con la función glm:

datos <- aplore3::chdage
attach(datos)
CHD <- as.integer(chd)-1  #Para codificar la variable *chd* como 0 y 1
modelo <- glm(CHD~age, family=binomial(link = "logit"), data=datos)
modelo$coefficients     # Todos los coeficientes
modelo$coefficients[1]  # Intercepto
modelo$coefficients[2]  # Pendiente

Con ambas funciones, se obtiene que: \(\hat{\delta}=\) -5.3094534 y \(\hat{\beta}=\) 0.1109211.

16.0.5 Solución parte (e)

Con la función summary del modelo glm obtenemos una salida en donde los valores de las estimaciones se pueden visualizar en la primera columna (llamada Estimate) de la lista llamada Coefficients (ver recuadro rojo en la figura 16.1):

datos <- aplore3::chdage
attach(datos)
CHD <- as.integer(chd)-1  #Para codificar la variable *chd* como 0 y 1
modelo <- glm(CHD~age, family=binomial(link = "logit"), data=datos)
summary(modelo)

Figure 16.1: Estimaciones de los parámetros logísticos. Fuente: Elaboración propia.

16.0.6 Solución parte (f)

Sabiendo que

\[\hat{p}_j = \hat{P}(\mbox{chd}=1 |\mbox{age}=x_j),\]

El modelo estimado se puede escribir de alguna de las dos maneras siguientes:

\[\mbox{Logit}(\hat{p}_j):= \ln\left(\frac{\hat{p}_j}{1-\hat{p}_j}\right) = -5.3094534 \;+\; 0.1109211 x_j, \qquad \qquad \hat{p}_j \;= \; \frac{e^{-5.3094534 \;+\; 0.1109211 x_j}} {1 + e^{-5.3094534 \;+\; 0.1109211 x_j}} \]

16.0.7 Solución parte (g)

La gráfica correspondiente se muestra en la figura 16.2. Se observa que hay una relación directa entre la edad y el riesgo de tener enfermedades coronarias:

edad <- seq(0, 100, 0.005)
delta <- modelo$coef[1]
beta <- modelo$coef[2]
logit <- delta+ beta*edad
p <- exp(logit)/(1+exp(logit))
D <- data.frame(edad,p)

ggplot(D, aes( x = edad)) +
geom_point(aes(y=p, colour=edad), size=1.5,color="red") +
labs(x="Edad", y="Riesgo de tener CHD", fill= "") + 
facet_wrap(. ~ "Riesgo de CHD versus Edad") +    
theme_bw(base_size = 12) 

Figure 16.2: Logit para los datos chdage

16.0.8 Solución parte (h)

La gráfica del logit estimado se puede visualizar en la Figura 16.3. Como es de esperar, hay una relación lineal y creciente entre el logit y la edad:

edad <- seq(0, 100, 0.005)
delta <- modelo$coef[1]
beta <- modelo$coef[2]
logit <- delta+ beta*edad
D <- data.frame(edad,logit)

ggplot(D, aes( x = edad)) +
geom_point(aes(y=logit, colour=edad), size=1.5,color="darkblue")  +
labs(x="Edad", y="Logit del riesgo de tener CHD", fill= "") + 
facet_wrap(. ~ "Logit del riesgo de CHD versus Edad") +    
theme_bw(base_size = 12) 

Figure 16.3: Logit para los datos chdage

16.0.9 Solución parte (i)

El logit estimado, para un sujeto con 50 años, es igual a 0.237. Observe que se puede obtener de varias maneras: reemplazando la edad en la ecuación (8.1) y con las funciones lsm() o glm():

\[\mbox{Logit}\Big(\hat{P}(\mbox{chd}=1 |\mbox{age}=50)\Big) \;= \; \hat{\delta} \;+\; 50 \hat{\beta} \; =\; -5.3094534 \;+\; (0.1109211)(50) \;= \; 0.2366\]

#1. Reemplazando:

edad <- 50
delta <- modelo$coef[1]
beta <- modelo$coef[2]
logit50 <- delta+ beta*edad

#2. Con la función lsm:

datos <- lsm::chdage
attach(datos)
modelo <- lsm(CHD~AGE, data=datos)
datos$logit <- modelo$Logit   #Incluye los logits en el data frame
logit50 <- subset(datos, age ==50)[1,4][1]

#3. Con la función glm:

datos <- aplore3::chdage
attach(datos)
CHD <- as.integer(chd)-1  #Para codificar la variable *chd* como 0 y 1
modelo <- glm(CHD~age, family=binomial(link = "logit"), data=datos)

#3a. Primera forma con glm:

datos$logit <- predict(modelo)   #Incluye los logits en el data frame 
logit50 <- subset(datos, age ==50)[1,5][1]

#3b. Segunda forma con glm:

datos$logit <- predict.glm(modelo)  
logit50 <- subset(datos, age ==50)[1,5][1]

#3c. Tercera forma con glm:

logit <- predict(modelo,newdata=data.frame(age=50)) #pj= 0.04347876 para age=20

16.0.10 Solución parte (j)

La proporción estimada de personas con presencia de enfermedades coronarias a la edad de 50 es igual a 0.559. Observe que se puede obtener de varias maneras: reemplazando la edad en la ecuación (9.1) y con las funciones lsm() o glm():

\[ \hat{P}(\mbox{chd}=1/\mbox{age}=50) \;= \; \frac{e^{\hat{\delta} \,+\, 50 \hat{\beta}}} {1 \;+\; e^{\hat{\delta} \,+\, 50\hat{\beta}}} \;= \; \frac{e^{-5.3094534 \;+\; (0.1109211)(50)}} {1 + e^{-5.3094534 \;+\; (0.1109211)(50)}} \;= \; \frac{e^{0.237}} {1 + e^{0.237}} \;= \; 0.5589\]

#1. Reemplazando:

edad <- 50
delta <- modelo$coef[1]
beta <- modelo$coef[2]
logit50 <- delta+ beta*edad
p50 <- exp(logit50)/(1+exp(logit50))

#2. Con la función lsm:

datos <- lsm::chdage
attach(datos)
modelo <- lsm(CHD~AGE, data=datos)
datos$logit <- modelo$Logit   #Incluye los logits en el data frame
logit50 <- subset(datos, age ==50)[1,4][1]
p50 <- exp(logit50)/(1+exp(logit50))

#3. Con la función glm:

datos <- aplore3::chdage
attach(datos)
CHD <- as.integer(chd)-1  #Para codificar la variable *chd* como 0 y 1
modelo <- glm(CHD~age, family=binomial(link = "logit"), data=datos)

#3a. Primera forma con glm:

datos$pj <- predict(modelo, type="response")   #Incluye los pj en el data frame 
p50 <- subset(datos, age ==50)[1,5][1]

#3b. Segunda forma con glm:

datos$pj <- predict.glm(modelo, type="response")  
p50 <- subset(datos, age ==50)[1,5][1]

#3c. Tercera forma con glm:

p50 <- predict(modelo,newdata=data.frame(age=50), type="response") 

16.0.11 Solución parte (k)

Los errores estándares de \(\hat{\delta}\) y \(\hat{\beta}\) son, respectivamente,

\[\hat{S}_{\hat{\delta}}= 1.13363, \qquad \hat{S}_{\hat{\beta}}= 0.0240593 \] Estas estimaciones se pueden obtener de varias maneras: con las funciones lsm() o glm::summary(). Con la función summary del modelo glm obtenemos una salida en donde los valores de las estimaciones de los errores estándares se pueden visualizar en la segunda columna (llamada Std. Error) de la lista llamada Coefficients (ver recuadro rojo en la figura 16.4):

#1. Con la función lsm:

datos <- lsm::chdage
attach(datos)
modelo <- lsm(CHD~AGE, data=datos)
StdError <- modelo$Std.Error

#2. Con la función glm::summary:

datos <- aplore3::chdage
attach(datos)
CHD <- as.integer(chd)-1  #Para codificar la variable *chd* como 0 y 1
modelo <- glm(CHD~age, family=binomial(link = "logit"), data=datos)
summary(modelo)

Figure 16.4: Estimaciones de los errores estándares. Fuente: Elaboración propia.

16.0.12 Solución parte (l)

La estimación del logaritmo de la función de máxima verosimilitud en el modelo logístico se puede obtener de varias maneras. Una es reemplazando los valores correspondientes en la ecuación (9.2):

\[{\cal L}(\hat{\alpha}) \;=\; \sum^{J}_{j=1}z_j\,\hat{g}_j \;-\; \sum^{J}_{j=1} n_j\ln \left[1+ e^{\hat{g}_j}\right] \;=\; -53.67655\]

En la expresión anterior, para un valor determinado \(x_j\) de la edad, \(\hat{g}_j\) es el logit estimado:

\[\hat{g}_j \;=\; \hat{\delta} \;+\; \hat{\beta} x_j = -5.3094534 \;+\; 0.1109211 \, x_j\]

datos <- lsm::chdage
attach(datos)
modelo <- lsm(CHD~AGE, data=datos)

#sum(modelo$z_j* modelo$Logit) - sum(modelo$n_j * ln(1+ exp(modelo$Logit)))
#modelo$Log_Lik_Logit
#modelo$n_j
dim(modelo$n_j)  #SALE NULL

## NULL

Las otras maneras para calcular \({\cal L}(\hat{\alpha})\) son con las funciones lsm() o glm::summary(). Con cualquier camino encontramos que \({\cal L}(\hat{\alpha})=\) -53.6765463.

Con la primera función se obtiene así:

#1. Con la función lsm:

datos <- lsm::chdage
attach(datos)
modelo <- lsm(CHD~AGE, data=datos)
modelo
modelo$Log_Lik_Logit

Con la función summary del modelo glm obtenemos una salida en donde se visualiza un valor (llamado Residual deviance) que se divide por -2 para obtener \({\cal L}(\hat{\alpha})\). En la parte (b) del ejemplo 2 de este mismo documento se explica este detalle. Este valor mencionado se puede observar en el último bloque de la salida (ver recuadro rojo en la figura 16.5):

#2. Con la función glm::summary:

datos <- aplore3::chdage
attach(datos)
CHD <- as.integer(chd)-1  #Para codificar la variable *chd* como 0 y 1
modelo <- glm(CHD~age, family=binomial(link = "logit"), data=datos)
summary(modelo)
modelo$deviance/-2

Figure 16.5: Estimación del Log de la función de verosimilitud en el modelo logístico. Fuente: Elaboración propia.

16.0.13 Solución parte (m)

La proporción estimada de personas con presencia de enfermedades coronarias a la edad de 50 es igual a 0.559. Observe que se puede obtener de varias maneras. Un camino es reemplazando la edad en la ecuación (9.1) :

\[ \hat{O}_{50} \;= \;\hat{O}(\mbox{cuando age=50}) \;= \; \frac{\hat{P}(\mbox{chd}=1/\mbox{age}=50)} {1 \;-\; \hat{P}(\mbox{chd}=1/\mbox{age}=50)} \;= \; \frac{0.5588765} {1 \;-\; 0.5588765 } \;= \; 1.266939\] Es decir, el odds estimado para individuos con edad de 50 es \(0,560/(1-0,560)=\) 1.266939. Esto significa que, cuando la persona tiene 50 años, la probabilidad de que tenga enfermedades coronarias es aproximadamente 1.27 veces la probabilidad de que no tenga.

Otra forma de hallar este valor estimado es utilizando las funciones lsm() o glm():

#1. Reemplazando:

edad <- 50
delta <- modelo$coef[1]
beta <- modelo$coef[2]
logit50 <- delta+ beta*edad
p50 <- exp(logit50)/(1+exp(logit50))
Odd50 <- p50/(1-p50)

#2. Con la función lsm:

datos <- lsm::chdage
attach(datos)
modelo <- lsm(CHD~AGE, data=datos)
datos$logit <- modelo$Logit   #Incluye los logits en el data frame
logit50 <- subset(datos, age ==50)[1,4][1]
p50 <- exp(logit50)/(1+exp(logit50))
Odd50 <- p50/(1-p50)

#3. Con la función glm:

datos <- aplore3::chdage
attach(datos)
CHD <- as.integer(chd)-1  #Para codificar la variable *chd* como 0 y 1
modelo <- glm(CHD~age, family=binomial(link = "logit"), data=datos)

#3a. Primera forma con glm:

datos$pj <- predict(modelo, type="response")   #Incluye los pj en el data frame 
p50 <- subset(datos, age ==50)[1,5][1]
Odd50 <- p50/(1-p50)

#3b. Segunda forma con glm:

datos$pj <- predict.glm(modelo, type="response")  
p50 <- subset(datos, age ==50)[1,5][1]
Odd50 <- p50/(1-p50)

#3c. Tercera forma con glm:

p50 <- predict(modelo,newdata=data.frame(age=50), type="response") 
Odd50 <- p50/(1-p50)

16.0.14 Solución parte (n)

La gráfica del odds estimado se puede visualizar en la Figura 16.6. Como es de esperar, hay una relación exponencial y creciente entre el odds y la edad:

edad <- seq(0, 100, 0.005)
delta <- modelo$coef[1]
beta <- modelo$coef[2]
logit <- delta+ beta*edad
p <- exp(logit)/(1+exp(logit))
odds <- p/(1-p)
D <- data.frame(edad,odds)

ggplot(D, aes( x = edad)) +
geom_point(aes(y=odds, colour=edad), size=1.5,color="darkgreen")  +
labs(x="Edad", y="Odds del riesgo de tener CHD", fill= "") + 
facet_wrap(. ~ "Odds del riesgo de CHD versus Edad") +    
theme_bw(base_size = 12) 

Figure 16.6: Logit para los datos chdage

16.0.15 Solución parte (o)

Por el teorema 12.1, cuando la edad se incrementa en 1 año, la razón odds estimada es 1.1173068:

\[\hat{OR} \;= \;e^{\hat{\beta}} \;=\; e^{0.11092} \;=\; 1.117307\]

Otra forma de hallar este valor estimado es utilizando las funciones lsm() o glm():

#2. Con la función lsm:

datos <- lsm::chdage
attach(datos)
modelo <- lsm(CHD~AGE, data=datos)

# 2a. Primera forma con lsm:

OR <- exp(modelo$coef[2])

#2b. Segunda forma con lsm:
OR <- modelo$ExpB[2,1]

#3. Con la función glm:

datos <- aplore3::chdage
attach(datos)
CHD <- as.integer(chd)-1  #Para codificar la variable *chd* como 0 y 1
modelo <- glm(CHD~age, family=binomial(link = "logit"), data=datos)
OR <- exp(modelo$coefficients[2])

17 Intervalos de confianza

En el documento Rpbus :: Modelos logísticos (intervalos de confianza) se ha desarrollado esta teoría, acompañado de las aplicaciones correspondientes.

18 Comparación de modelos

Esta sección se refiere a resaltar que hay estadísticos para distintas pruebas de comparación de modelos:

1. \(H_0\): Logístico versus \(H_1\): Completo.

2. \(H_0\): Nulo versus \(H_1\): Logístico.

3. \(H_0\): Logístico versus \(H_1\): Saturado.

4. \(H_0\): Submodelo versus \(H_1\): Logístico.

5. \(H_0\): Submodelo con una variable explicativa menos versus \(H_1\): Logístico.

En el documento Rpbus :: Modelos logísticos (pruebas de hipótesis) se ha desarrollado esta teoría, acompañada con las aplicaciones correspondientes.

19 Matriz de confusión

En el campo del Aprendizaje automático (Machine Learning), una matriz de confusión (también conocida como matriz de error), es una tabla de contingencia, con dos dimensiones (“real, observado” y “predicho”) que permite la visualización del rendimiento de un algoritmo, típicamente uno de aprendizaje supervisado. Por lo general, en aprendizaje no supervisado se denomina matriz de emparejamiento (en inglés: matching matrix). Cada fila de la matriz representa las ocurrencias en la clase observada, mientras que cada columna representa las ocurrencias en la clase predicha, o viceversa. El nombre proviene del hecho de que facilita determinar si el sistema está confundiendo dos clases (es decir, comúnmente etiquetado incorrectamente una como otra). En la Tabla 19.1 se muestra la estructura general de una matrix de confusión:

Table 19.1: Matrix de confusión.

1	2	3	4	5	6
Condiciones	Positiva predicha (\(y=1\))	Negativa predicha (\(y=0\))	Total	Estadísticos	Estadísticos
Positiva observada (\(y=1\))	Verdadero positivo (\(a\))	Falso negativo (\(c\))	\(a+c\)	Sensitividad = \(a/(a+c)\)	Prevalencia \(=(a+c)/n\)
Negativa observada (\(y=0\))	Falso positivo (\(b\))	Verdadero negativo (\(d\))	\(b+d\)	Especificidad = \(d/(b+d)\)	1- Especificidad \(=b/(b+d)\)
Total	\(a+b\)	\(c+d\)	Tamaño muestral (\(n=a+b+c+d)\)
Estadísticos	Valor Predictivo positivo \(=a/(a+b)\)	Valor predictivo negativo \(=d/(c+d)\)		Exactitud = \((a+d)/n\)

Teniendo en cuenta la tabla anterior, se pueden definir los siguientes conceptos:

1. La exactitud ( en inglés: accuracy) es la proporción de clasificaciones correctas.

\[Exactitud \; =\; \frac{a+d}{n}\]

2. La sensitividad o tasa de verdaderos positivos (en inglés: sensivity, true positive rate, recall or hit rate) es la proporción de casos positivos observados que son correctamente clasificados:

\[ Sensitividad \;=\; \frac{a}{a+c} \]

3. La especificidad o tasa de verdaderos negativos (en inglés: specifity, selectivity or true negative rate) es la proporción de casos negativos observados que son correctamente clasificados:

\[ Especificidad \;=\; \frac{d}{b+d} \]

4. La tasa de falsos positivos (1- especificidad) es la proporción de casos negativos observados que son incorrectamente clasificados:

\[ 1- Especificidad \;=\; \frac{b}{b+d} \]

5. La precisión o valor predictivo positivo VPP (en inglés: precision, positive predictive value) es la proporción de casos positivos que son correctamente clasificados:

\[ VPP \;=\; \frac{a}{a+b} \]

6. El valor predictivo negativo VPN (en inglés: negative predictive value) es la proporción de casos negativos que son correctamente clasificados:

\[ VPN \;=\; \frac{d}{c+d} \]

5. La prevalencia (en inglés, prevalence) es la proporción de casos positivos observados:

\[ Prevalencia \;=\; \frac{a+c}{n} \]

6. La exactitud balanceada (en inglés: balanced accuracy) se define como:

\[\mbox{Exactitud balanceada} \;=\; \frac{\mbox{sensitividad} + \mbox{especificidad}}{2}\]

7. La tasa de detección (en inglés: Detection rate) se define como la proporción de casos positivos observados y predichos:

\[\mbox{Tasa de detección} \;=\; \frac{a}{n}\] 8. La prevalencia de detección (en inglés: Detection prevalence) se define como la proporción de casos positivos predichos:

\[\mbox{Prevalencia de detección} \;=\; \frac{a+b}{n}\]

20 Indices de concordancias

Existen varios índices de concordancia propuestos y toman valores entre 0 (total desacuerdo) y 1 (máximo acuerdo). Un ejemplo es la exactitud (definido anteriormente), pero el más usado es el denominado índice kappa de Cohen (\(\kappa\)) y se define como:

\[\kappa \;=\; \frac{Exactitud - P_E}{1- P_E}\]

donde \(P_E\) es la proporción de acuerdos esperados. Está dado por:

\[\begin{eqnarray*} P_E &=& \frac{(\mbox{Total fila y=1)(Total columna y=1)} \;+\;\mbox{(Total fila y=0)(Total columna y=0)}}{\mbox{(Tamaño de la muestra)}^2} \\ && \\ &=& \frac{(a+c)(a+b) + (b+d)(c+d)}{n^2} \end{eqnarray*}\]

Se ha utilizado la siguiente escala de valoración del grado de acuerdo para \(\kappa\):

1. \(\kappa \leq 0.00\): sin acuerdo.

2. \(0.00 < \kappa \leq 0.20\): insignificante.

3. \(0.20 <\kappa \leq 0.40\): discreto.

4. \(0.40 < \kappa \leq 0.60\): moderado.

5. \(0.60 < \kappa \leq 0.80\): sustancial.

6. \(0.80 < \kappa \leq 1.00\): casi perfecto.

21 Curva ROC

Una curva característica de funcionamiento del receptor, o curva ROC (en inglés: receiver operating characteristic curve), es un gráfico (como la que se muestra en la Figura 21.1a) que ilustra la capacidad de diagnóstico de un sistema clasificador binario a medida que varía su umbral de discriminación. La curva ROC se construye sobre un espacio (llamado espacio ROC) trazando la tasa de verdaderos positivos (sentitividad) frente a la tasa de falsos positivos (1-especificidad) para varios valores de ciertos umbrales (véase también la Figura 21.1b).

Figure 21.1: Ejemplo de una curva y de un espacio ROC.

Cada resultado de predicción de una matriz de confusión representa un punto en el espacio ROC. El mejor método de predicción posible produciría un punto en la esquina superior izquierda o coordenada (0,1) del espacio ROC, que representa el 100% de sensibilidad (sin falsos negativos) y el 100% de especificidad (sin falsos positivos). El punto (0,1) también se denomina clasificación perfecta. Una suposición aleatoria daría un punto a lo largo de una línea diagonal (la llamada línea de no discriminación) desde la esquina inferior izquierda hasta la esquina superior derecha (independientemente de las tasas base positivas y negativas). Es importante resaltar que la diagonal divide el espacio ROC. Los puntos por encima de la diagonal representan buenos resultados de clasificación y los puntos debajo de la línea representan malos resultados.

El área bajo la curva ROC (AUC, en inglés: Area Under the Curve) provee una medida de la habilidad del modelo para discriminar entre aquellos sujetos que experimentan la variable de respuesta de interés versus aquéllos que no lo hacen. Como regla general:

1. Si \(AUC=0.5\), entonces, no hay discriminación.

2. Si \(0.7 \leq AUC < 0.8\), la discriminación se considera aceptable.

3. Si \(0.8 \leq AUC < 0.9\), la discriminación se considera excelente.

4. Si \(AUC \geq 0.9\), la discriminación se considera sobresaliente (en la práctica, es extremadamente inusual observar áreas con esta condición).

22 Paquetes R para curvas ROC

En la CRAN hemos encontrado 93 paquetes que contienen el término ROC. De ellos, hay 52 que tienen un puntaje (score) mayor que 250. En la Tabla 22.1 se muestran los 6 primeros.

library(pkgsearch)
library(dplyr)

Paq_ROC <-  pkg_search(query="ROC",size=200)
nrow(Paq_ROC)

Paq_ROC_Red <- Paq_ROC %>% 
               filter(maintainer_name != "ORPHANED", score > 250) %>%
               select(score, package, downloads_last_month) %>%
               arrange(desc(downloads_last_month))
head(Paq_ROC_Red)
nrow(Paq_ROC_Red)

Table 22.1: Paquetes CRAN relacionados con ROC.

1	2	3
score	package	downloads_last_month
4512.488	caTools	133410
10160.324	pROC	83818
952.7672	ROCR	63403
2401.1794	PRROC	9299
250.35356	riskRegression	7305
1848.5184	cvAUC	4448

La figura 22.1 utiliza el paquete dlstats con el fin de mostrar el historial de descargas de los seis paquetes que se seleccionaron en la tabla anterior.

library(ggplot2)
library(dlstats)
Lista <- c("caTools", "pROC", "ROCR","PRROC", "riskRegression", "sROC")
PaquetesROC <- cran_stats(Lista)

ggplot(PaquetesROC, aes(end, downloads, group=package, color=package)) +
  geom_line() + geom_point(aes(shape=package)) +
  scale_y_continuous(trans = 'log2')

Figure 22.1: Historial de los 6 paquetes tops de la CRAN que utilizan ROC

23 Indice de Youden

Las pruebas de diagnósticos que producen rangos o puntuaciones pueden ser analizadas con un punto de corte o umbral de decisión (\(c\)), obteniendo así un clasificador binario: si el marcador está por encima o es igual al punto de corte produce un Si (clase positiva) y, si no es así, produce un No (clase negativa). Además, la sensibilidad y la especificidad de una prueba diagnóstica son usualmente utilizadas, de forma simultánea, como una medida conjunta del comportamiento de esa prueba. Esto se debe a que son complementarias. En general, al disminuir \(c\), la fracción de verdaderos positivos aumenta, dando lugar a que la fracción de verdaderos negativos disminuya . En esta situación, hay que alcanzar un compromiso aceptable. Una de las soluciones propuestas es seleccionar el punto de corte que maximice la diferencia entre las fracciones de verdaderos positivos y falsos positivos. El valor máximo de esta cantidad es el llamado índice de Youden (\(J\)) y el punto de corte \(c\) en el punto de la curva ROC correspondiente a este índice es, a menudo, seleccionado como el punto de corte óptimo de la prueba. Este índice está definido por \[J \;= \max(\mbox{Sensitividad} \;+\; \mbox{Especificidad} \;-\; 1)\]

El índice fue sugerido por W.J. Youden en 1950 como una forma de resumir el desempeño de una prueba de diagnóstico. Su valor varía de 0 a 1 (inclusive):

1. Es cero cuando una prueba de diagnóstico da la misma proporción de resultados positivos para grupos con y sin la enfermedad. Es decir, cuando la prueba es inútil.

2. Un valor de 1 indica que no hay falsos positivos ni falsos negativos. Es decir, la prueba es perfecta.

24 Ejemplo 2: enunciado

Considere los datos del archivo chdage. Suponga que se quiere analizar un modelo de regresión logística, considerando a chd como variable dependiente y age como independiente.

1. Calcule tanto las probabilidades predichas como los valores predichos de la variable de respuesta y reúnalos en el data frame chdage.

2. Teniendo en cuenta el data frame de (a), construya la matriz de confusión correspondiente.

3. Determine la exactitud del modelo.

4. Aplique la librería caret para hallar la matriz de confusión y los estadísticos relacionados (exactitud, sensitividad, especificidad, etc.).

5. Aplique la librería gmodels para hallar la matriz de confusión y los estadísticos relacionados (exactitud, sensitividad, especificidad, etc.).

6. Construya la curva ROC y halle el AUC.

7. Construya un data frame que contenga las sensitividades, las especificidades y los umbrales correspondientes.

8. Encuentre el índice de Youden.

9. Construya la curva ROC y halle el AUC, pero cuando la tasa de falsos positivos se encuentra entre 10% y 30%.

25 Ejemplo 2: Solución

25.0.1 Solución parte (a)

Abajo se muestra una parte del data frame solicitado:

library(lsm)
datos <- lsm::chdage
attach(datos)
modelo <- lsm(CHD~AGE, data=datos)


#2. Con la función glm:

datos <- aplore3::chdage
attach(datos)
CHD <- as.integer(chd)-1  #Para codificar la variable *chd* como 0 y 1
modelo <- glm(CHD~age, family=binomial(link = "logit"), data=datos)

datos$Prob_pred <- round(predict(modelo, newdata = datos, "response"),4)

datos$CHD_pred <- ifelse(datos$Prob_pred > 0.5, "Yes", "No")

y_pred <- as.factor(datos$CHD_pred)
y_obs <- as.factor(datos$chd)
datos

Table 25.1: Data frame con los predichos.

id	age	agegrp	chd	Prob_pred	CHD_pred
1	20	20-39	No	0.0435	No
2	23	20-39	No	0.0596	No
3	24	20-39	No	0.0662	No
4	25	20-39	No	0.0733	No
5	25	20-39	Yes	0.0733	No
:	:			:	:
:	:			:	:
96	63	60-69	Yes	0.8427	Yes
97	64	60-69	No	0.8569	Yes
98	64	60-69	Yes	0.8569	Yes
99	65	60-69	Yes	0.8699	Yes
100	69	60-69	Yes	0.9125	Yes

25.0.2 Solución parte (b)

La matriz de confusión correspondiente es:

Tabla <-  table(y_obs, y_pred)
Tabla

Table 25.2: Matriz de confusión.

	1	2	3
	y_pred=No	y_pred=Yes	Total
y_obs=No	45	12	57
y_obs=Yes	14	29	43
Total	59	41	100

25.0.3 Solución parte (c)

La exactitud es del 74% y se puede calcular de dos maneras:

# Con el data frame de la parte (a):
E <- mean(y_pred == y_obs)*100

# Con la matriz de confusión de la parte (b):
E <- round((sum(diag(Tabla))/sum(Tabla))*100,2)

25.0.4 Solución parte (d)

Con la función confusionMatrix de la librería caret se pueden hallar la matriz de confusión y los estadísticos solicitados (véase la figura 25.1):

library(caret)
confusionMatrix(data=y_pred, reference = y_obs,  positive="Yes")

Figure 25.1: Matriz de confusión con el paquete caret. Fuente: Elaboración propia.

25.0.5 Solución parte (e)

Con la función CrossTable de la librería gmodels se pueden hallar la matriz de confusión y también los estadísticos solicitados (véase la figura 25.2):

library(gmodels)
CrossTable(y_pred, y_obs)

Figure 25.2: Matriz de confusión con el paquete gmodels. Fuente: Elaboración propia.

25.0.6 Solución parte (f)

La curva ROC solicitada es:

datos <- aplore3::chdage
attach(datos)
CHD <- as.integer(chd)-1  #Para codificar la variable *chd* como 0 y 1
modelo <- glm(CHD~age, family=binomial(link = "logit"), data=datos)

y_obs <- as.factor(datos$chd)

library("pROC")

par(pty = "s") # Hacer cuadrado el espacio ROC

roc_graph <- roc(y_obs, modelo$fitted.values, plot = TRUE, legacy.axes = TRUE,
    percent = TRUE, ylab = "Sensitividad \n (Porcentaje de verdaderos positivos)",
    xlab = "1- especificidad \n (Porcentaje de falsos positivos)", col = "darkblue", lwd = 2,
    print.auc = TRUE, auc.polygon = TRUE, auc.polygon.col = "aliceblue", grid=TRUE)
coordenada <- pROC::coords(roc_graph,"best",ret=c("threshold","specificity","sensitivity"))
points(x=coordenada$specificity, y=coordenada$sensitivity, pch=0,col="red")
roc_graph

## 
## Call:
## roc.default(response = y_obs, predictor = modelo$fitted.values,     percent = TRUE, plot = TRUE, legacy.axes = TRUE, ylab = "Sensitividad \n (Porcentaje de verdaderos positivos)",     xlab = "1- especificidad \n (Porcentaje de falsos positivos)",     col = "darkblue", lwd = 2, print.auc = TRUE, auc.polygon = TRUE,     auc.polygon.col = "aliceblue", grid = TRUE)
## 
## Data: modelo$fitted.values in 57 controls (y_obs No) < 43 cases (y_obs Yes).
## Area under the curve: 79.99%

25.0.7 Solución parte (g)

El data frame solicitado se muestra abajo:

roc_datos <- data.frame(Sensitividad =roc_graph$sensitivities/100,
                     Especificidad = roc_graph$specificities/100,
                     Umbrales = roc_graph$thresholds)
head(roc_datos)

Table 25.3: Data frame con las sensitividades, las especificidades y los umbrales.

Posicion	Sensitividad	Especificidad	Umbrales
1	1	0	-Inf
2	1	0.018	0.052
3	1	0.035	0.063
4	1	0.053	0.07
5	0.977	0.07	0.077
:	:	:	:
21	0.791	0.684	0.381
:	:	:	:
40	0.093	0.982	0.835
41	0.07	0.982	0.85
42	0.047	1	0.863
43	0.023	1	0.891
44	0	1	Inf

25.0.8 Solución parte (h)

El índice de Youden se encuentra en la posición 21 del data frame construido en el inciso (h) y es 0.381.

#Primera forma:
roc_datos[which.max(roc_datos$Sensitividad+roc_datos$Especificidad-1),]

#Segunda forma:
pROC::coords(roc_graph, x="best", input="threshold", best.method="youden")

#Tercera forma:
pROC::coords(roc_graph,"best",ret=c("threshold","specificity","sensitivity"))

Figure 25.3: Indice de Youden. Fuente: Elaboración propia.

25.0.9 Solución parte (i)

Cuando nos interesa una parte de la curva ROC que solo permita pequeños numeros de falsos positivos (por ejemplo, cuando la tasa de falsos positivos se encuentra entre 10% y 30%), procedemos de la siguiente manera:

datos <- aplore3::chdage
attach(datos)
CHD <- as.integer(chd)-1  #Para codificar la variable *chd* como 0 y 1
modelo <- glm(CHD~age, family=binomial(link = "logit"), data=datos)

y_obs <- as.factor(datos$chd)

library("pROC")
#library("randomForest")

par(pty = "s") # Hacer cuadrado el espacio ROC

roc_graph <- roc(y_obs, modelo$fitted.values, plot = TRUE, legacy.axes = TRUE,
    percent = TRUE, ylab = "Sensitividad \n (Porcentaje de verdaderos positivos)",
    xlab = "1- especificidad \n (Porcentaje de falsos positivos)", col = "darkblue", lwd = 2,
    print.auc = TRUE, print.auc.x =57, partial.auc = c(90, 70), # en términos de especificidad
    auc.polygon = TRUE, auc.polygon.col = "aliceblue", grid=TRUE)
coordenada <- pROC::coords(roc_graph,"best",ret=c("threshold","specificity","sensitivity"))
points(x=coordenada$specificity, y=coordenada$sensitivity, pch=0,col="red")
roc_graph

## 
## Call:
## roc.default(response = y_obs, predictor = modelo$fitted.values,     percent = TRUE, plot = TRUE, legacy.axes = TRUE, ylab = "Sensitividad \n (Porcentaje de verdaderos positivos)",     xlab = "1- especificidad \n (Porcentaje de falsos positivos)",     col = "darkblue", lwd = 2, print.auc = TRUE, print.auc.x = 57,     partial.auc = c(90, 70), auc.polygon = TRUE, auc.polygon.col = "aliceblue",     grid = TRUE)
## 
## Data: modelo$fitted.values in 57 controls (y_obs No) < 43 cases (y_obs Yes).
## Partial area under the curve (specificity 90%-70%): 13.02%

26 Ejercicios

Para la solución de los siguientes ejercicios, téngase en cuenta los siguientes comentarios:

• Todos los datos mencionados aparecen en los paquetes mencionados en este documento.

• Siempre debe detallar el análisis del conjunto de datos (con las variables especificadas) basado en lo explicado en este documento.

• Verifique cómo se obtienen las estimaciones correspondientes, los logaritmos de las funciones de máxima verosimilitud, estimaciones e intervalos de confianza para \(p_j\), ODDS, razones ODDS, intercepto, pendiente, etc.

26.0.1 Ejercicios 1 a 3

1. Demuestre los teoremas: (a) 5.1; (b) 6.2; (c) 7.1; (d) 7.2; (e) 9.1

2. Haga un listado de los paquetes de R que estimen el logaritmo de la función de máxima verosimilitud en los modelos completo y nulo.

3. Haga un listado de los paquetes de R que, en el caso binario, estimen el logaritmo de la función de máxima verosimilitud en los modelos saturado y logístico.

26.0.2 Ejercicio 4

4. Los datos ICU corresponden a una muestra de 200 sujetos que hicieron parte de un estudio de supervivencia de pacientes que fueron remitidos a una unidad de cuidados intensivos (intensive care unit - ICU). La meta principal de este estudio fue desarrollar un modelo de regresión logística para predecir la probabilidad de supervivencia de estos pacientes en el hospital y estudiar los factores de riesgos asociados con el índice de mortalidad ICU. En estos datos tome a la variable AGE como independiente y STA como dependiente.

1. Escriba la ecuación general para el modelo de regresión logístico de STA contra AGE y para el logit transformado de este modelo. ¿Qué características de STA nos pone a pensar que debamos considerar el modelo de regresión logística en vez del usual modelo de regresión lineal para describir la relación entre STA y AGE?

2. Forme un diagrama de dispersión de STA contra AGE.

3. Usando los intervalos [15,24], [25,34], [35,44], [45,54], [55,64], [65,74], [75,84], [85,94] para AGE, calcule la media de STA de los sujetos dentro de cada intervalo. Grafique estos valores de la media de STA contra el punto medio del intervalo de AGE usando el mismo conjunto de ejes que se utilizaron en la parte (b).

4. Escriba una expresión para la función de verosimilitud y del logaritmo de esta función para el modelo de regresión logístico de (a) usando los 200 datos no agrupados. Obtenga una expresión para las dos ecuaciones de verosimilitud.

5. Obtenga las estimaciones de los parámetros del modelo de regresión logístico de (a). Usando estas estimaciones, escriba las correspondientes ecuaciones para los valores ajustados. Grafique la ecuación para los valores ajustados utilizando los mismos ejes como en (b) y (c).

6. Resuma (describa en palabras) los resultados presentados en la gráfica obtenida en (b), (c) y (e).

7. Usando los resultados de la parte (e), verifique la significancia del coeficiente de AGE. ¿Qué supuestos se necesitan para realizar dicha prueba?

8. Usando los resultados de (e), halle un intervalo del 95% de confianza para la pendiente y la constante. Escriba una interpretación con respecto al intervalo encontrado para la pendiente.

9. Obtenga la matriz de covarianzas estimada para el modelo en (e). Calcule el logit y la probabilidad logística estimada para una persona de 60 años. Calcule un intervalo del 95% de confianza para el logit y la probabilidad logística estimada. Interprete la probabilidad estimada y su intervalo de confianza.

10. Obtenga el logit estimado y su error estándar para cada persona en el estudio ICU. Grafique el logit estimado y los límites del intervalo del 95 % de confianza versus AGE para cada persona. Explique (en palabras) similaridades y diferencias entre las apariencias de esta gráfica y una gráfica de una gráfica de un modelo de regresión ajustado y sus límites del intervalo del 95 % de confianza.

11. Realize las siguientes pruebas de comparación de modelos resumiendo en una tabla las pruebas realizadas, el valor y la distribución muestral del estadístico de prueba, los grados de libertad, el P-valor y su decisión.

• Nulo vs Logístico.

• Logístico vs Completo.

• Logístico vs Saturado

26.0.3 Ejercicios 5 a 7

5. Considere los datos ICU. Repita el ejercicio 4 utilizando la variable TYP (como variable dependiente) en vez de STA.

6. Considere los datos ICU. Repita todos los análisis realizados en este documento, pero considerando ahora las variables AGE (como variable independiente) y STA (como variable dependiente).

7. Considere los datos ICU. Haga el análisis correspondiente tomando a STA como variable dependiente y a AGE, SYS y HRA como independientes.

26.0.4 Ejercicios 8 a 9

8. Los datos UIS se recogieron con el propósito de comparar dos programas de tratamiento A y B para reducir el abuso de la droga y prevenir sus riesgos. La descipción de los datos se puede ver también aquí. Detalle el análisis para estos datos, tomando a DFREE como variable dependiente y AGE, BECK y NDRUGTX como variables independientes.

9. Los datos PROS corresponden a un estudio realizado pacientes con cáncer de próstata para determinar si las variables medidas en un examen básico pueden ser usadas para predecir si el tumor ha penetrado la cápsula prostática. Los datos fueron recogidos teniendo en cuenta 380 individuos, 153 de los cuales tuvieron un cáncer que penetró la cápsula prostática. En estos datos, una variable que se pensó que era particularmente predictiva para la penetración de cápsula es el nivel de antígeno prostático, PSA. Repita los pasos del ejercicio 4 usando CAPSULE como variable dependiente y utilize para PSA, los intervalos [0.0; 2.4], [2.5; 4.4], [4.5; 6.4], [6.5; 8.4], [8.5; 10.4], [10.5; 12.4], [12.5; 20.4], [20.5; 140].

26.0.5 Ejercicio 10

10. De todas las variables que aparecen en los datos PROS sólo considere a CAPSULE (como variable dependiente) y PSA (como variable independiente).

1. Responda:

• ¿Cuál es la ecuación para el modelo de regresión logística?

• ¿Cuál es la ecuación para riesgo estimado por este modelo?

• ¿Qué características de la variable dependiente nos conduce a considerar la regresión logística como más apropiada que el modelo de regresión lineal para describir la relación entre las dos variables mencionadas anteriormente?

2. Calcule:

• \(\mathcal{L}(\widetilde{p})\) en el modelo saturado.

• \(\mathcal{L}(\widehat{\alpha})\) en el modelo logístico.

3. Halle las estimaciones de los siguientes parámetros e interprételos (justifique en forma clara y precisa todas sus afirmaciones):

• \(P(CAPSULE=1 \, / \, PSA=11.2\) mg/ml\()\).

• \(P(CAPSULE=0 \, / \, PSA=11.2\) mg/ml\()\).

• El odds cuando PSA=11.2.

• La razón odds OR.

4. Realize las siguientes pruebas de comparación de modelos resumiendo en una tabla las pruebas realizadas, el valor y la distribución muestral del estadístico de prueba, los grados de libertad, el P-valor y su decisión.

• Nulo vs Logístico.

• Logístico vs Completo.

• Logístico vs Saturado

26.0.6 Ejercicios 11 a 13

11. Considere los datos PROS. Realize un análisis de regresión logística tomando a CAPSULE como variable dependiente y VOL como variable independiente.

12. Considere los datos PROS. Realize un análisis de regresión logística tomando a CAPSULE como variable dependiente y AGE como variable independiente.

13. Considere los datos PROS, tomando a CAPSULE como variable dependiente y AGE, PSA y VOL como variables independientes.

26.0.7 Ejercicio 14

14. Los datos LOWBWT corresponden a un estudio realizado para identificar factores de riesgos asociados a nacimientos de bebés con bajo peso (peso menor que 2.500 gramos). Los datos fueron recogidos teniendo en cuenta 189 mujeres, 59 de las cuales tuvieron bebés con bajo peso y 130 de las cuales tuvieron bebés con peso normal. De todas las variables que aparecen sólo considere a LOW (como variable dependiente) y LWT (como variable independiente).

1. Responda:

• ¿Cuál es la ecuación para el modelo de regresión logística?

• ¿Cuál es la ecuación para el riesgo estimado por este modelo?

• ¿Qué características de la variable dependiente nos conduce a considerar la regresión logística como más apropiada que el modelo de regresión lineal para describir la relación entre las dos variables mencionadas anteriormente?

2. Calcule:

• \(\mathcal{L}(\widetilde{p})\) en el modelo saturado.

• \(\mathcal{L}(\widehat{\alpha})\) en el modelo logístico.

3. Halle las estimaciones para los siguientes parámetros e interprételos (justifique en forma clara y precisa todas sus afirmaciones):

• \(P(LOW=1 \, / \, LWT=100.3\) libras\()\).

• \(P(LOW=0 \, / \, LWT=100.3\) libras\()\).

• El odds cuando LWT=100.3.

• La razón odds OR.

4. Realize las siguientes pruebas de comparación de modelos resumiendo en una tabla las pruebas realizadas, el valor y la distribución muestral del estadístico de prueba, los grados de libertad, el P-valor y su decisión.

• Nulo vs Logístico.

• Logístico vs Completo.

• Logístico vs Saturado

26.0.8 Ejercicios 15 a 18

15. Considere los datos LOWBWT, tomando a LOW como variable dependiente y AGE como variable independiente.

16. Considere los datos LOWBWT, tomando a LOW como variable dependiente y LWT como variable independiente.

17. Considere los datos LOWBWT, tomando a LOW como variable dependiente y BWT como variable independiente.

18. Considere los datos LOWBWT, tomando a LOW como variable dependiente y AGE, LWT y BWT como variables independientes.

26.0.9 Ejercicios 19 a 21

19. Considere los datos LOWBWTM11, tomando a LOW como variable dependiente y AGE como variable independiente.

20. Considere los datos LOWBWTM11, tomando a LOW como variable dependiente y LWT como variable independiente.

21. Considere los datos LOWBWTM11, tomando a LOW como variable dependiente y AGE y LWT como variables independientes.

26.0.10 Ejercicios 22

22. En los datos CLSLOWBWT, entre otras variables, una que los físicos consideraron importante para el control del peso del bebé (variable dependiente LOW) fue el peso de la madre durante el último periodo menstrual (LWT). Repita los pasos del ejercicio 4, pero para la parte (c) utilize los intervalos: [80,99], [100,109], [110,114], [115,119], [120,124], [125,129], [130,250]. La gráfica en la parte (c) no parece en forma de \(S\). La razón principal es que el rango de los valores graficados está aproximadamente entre 0.2 y 0.56. Explique por qué un modelo para la probabilidad de LOW como una función de LWT pudiese ser el modelo de regresión logística.

26.0.11 Ejercicios 23 a 26

23. Considere los datos CLSLOWBWT, tomando a LOW como variable dependiente y AGE como variable independiente.

24. Considere los datos CLSLOWBWT, tomando a LOW como variable dependiente y LWT como variable independiente.

25. Considere los datos CLSLOWBWT, tomando a LOW como variable dependiente y BWT como variable independiente.

26. Considere los datos CLSLOWBWT, tomando a LOW como variable dependiente y AGE, LWT y BWT como variables independientes.

Bibliografía

Consultar el documento RPubs :: Regresión logística (bibliografía).

 

 

If you found any ERRORS or have SUGGESTIONS, please report them to my email. Thanks.