1. Introducción

Comenzaremos examinando algunos resúmenes numéricos y gráficos de los datos del mercado, que es parte de la librería ISLR. Este conjunto de datos consiste en retornos porcentuales para el índice bursátil S&P 500 durante 1,250 días, desde el comienzo de 2001 hasta finales de 2005. Para cada fecha, hemos registrado los retornos porcentuales para cada uno de los cinco días de negociación anteriores, Lag1 a Lag5. También hemos registrado el Volume (el número de acciones negociadas en el día anterior, en miles de millones), Today (el rendimiento porcentual en la fecha en cuestión) y Direction (si el mercado estaba Arriba o Abajo en esta fecha)

library(ISLR)
names(Smarket)

## [1] "Year"      "Lag1"      "Lag2"      "Lag3"      "Lag4"      "Lag5"     
## [7] "Volume"    "Today"     "Direction"

dim(Smarket)

## [1] 1250    9

summary(Smarket)

##       Year           Lag1                Lag2                Lag3          
##  Min.   :2001   Min.   :-4.922000   Min.   :-4.922000   Min.   :-4.922000  
##  1st Qu.:2002   1st Qu.:-0.639500   1st Qu.:-0.639500   1st Qu.:-0.640000  
##  Median :2003   Median : 0.039000   Median : 0.039000   Median : 0.038500  
##  Mean   :2003   Mean   : 0.003834   Mean   : 0.003919   Mean   : 0.001716  
##  3rd Qu.:2004   3rd Qu.: 0.596750   3rd Qu.: 0.596750   3rd Qu.: 0.596750  
##  Max.   :2005   Max.   : 5.733000   Max.   : 5.733000   Max.   : 5.733000  
##       Lag4                Lag5              Volume           Today          
##  Min.   :-4.922000   Min.   :-4.92200   Min.   :0.3561   Min.   :-4.922000  
##  1st Qu.:-0.640000   1st Qu.:-0.64000   1st Qu.:1.2574   1st Qu.:-0.639500  
##  Median : 0.038500   Median : 0.03850   Median :1.4229   Median : 0.038500  
##  Mean   : 0.001636   Mean   : 0.00561   Mean   :1.4783   Mean   : 0.003138  
##  3rd Qu.: 0.596750   3rd Qu.: 0.59700   3rd Qu.:1.6417   3rd Qu.: 0.596750  
##  Max.   : 5.733000   Max.   : 5.73300   Max.   :3.1525   Max.   : 5.733000  
##  Direction 
##  Down:602  
##  Up  :648  
##            
##            
##            
##

La función cor() produce una matriz que contiene todas las correlaciones por pares entre los predictores en un conjunto de datos. Tenemos que quitar la variable Dirección del análisis ya que es cualitativa.

cor(Smarket[,-9])

##              Year         Lag1         Lag2         Lag3         Lag4
## Year   1.00000000  0.029699649  0.030596422  0.033194581  0.035688718
## Lag1   0.02969965  1.000000000 -0.026294328 -0.010803402 -0.002985911
## Lag2   0.03059642 -0.026294328  1.000000000 -0.025896670 -0.010853533
## Lag3   0.03319458 -0.010803402 -0.025896670  1.000000000 -0.024051036
## Lag4   0.03568872 -0.002985911 -0.010853533 -0.024051036  1.000000000
## Lag5   0.02978799 -0.005674606 -0.003557949 -0.018808338 -0.027083641
## Volume 0.53900647  0.040909908 -0.043383215 -0.041823686 -0.048414246
## Today  0.03009523 -0.026155045 -0.010250033 -0.002447647 -0.006899527
##                Lag5      Volume        Today
## Year    0.029787995  0.53900647  0.030095229
## Lag1   -0.005674606  0.04090991 -0.026155045
## Lag2   -0.003557949 -0.04338321 -0.010250033
## Lag3   -0.018808338 -0.04182369 -0.002447647
## Lag4   -0.027083641 -0.04841425 -0.006899527
## Lag5    1.000000000 -0.02200231 -0.034860083
## Volume -0.022002315  1.00000000  0.014591823
## Today  -0.034860083  0.01459182  1.000000000

Como era de esperar, las correlaciones entre las variables de retraso y los rendimientos actuales son cercanas a cero. En otras palabras, parece haber poca correlación entre los retornos de hoy y los retornos de días anteriores. La única correlación sustancial es entre Year y Volumen. Al trazar los datos se ve que el volumen aumenta con el tiempo. En otras palabras, el número promedio de acciones negociadas diariamente aumentó de 2001 a 2005.

attach(Smarket)
plot(Volume)

2. Regresión Logística

A continuación, ajustaremos un modelo de regresión logística para predecir Direction usando Lag1 a Lag5 y Volumen. La función glm() se ajusta a modelos lineales glm () generalizados, una clase de modelos que incluye regresión logística. El modelo lineal generalizado de sintaxis de la función glm() es similar al de lm(), excepto que debemos incluirel argumento family=binomial para decirle a R que ejecute una regresión logística en lugar de otro tipo de modelo lineal generalizado.

contrasts(Smarket$Direction)

##      Up
## Down  0
## Up    1

glm.fit<-glm(Direction~Lag1+Lag2+Lag3+Lag4+Lag5+Volume,data=Smarket,family=binomial)
summary(glm.fit)

## 
## Call:
## glm(formula = Direction ~ Lag1 + Lag2 + Lag3 + Lag4 + Lag5 + 
##     Volume, family = binomial, data = Smarket)
## 
## Coefficients:
##              Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.126000   0.240736  -0.523    0.601
## Lag1        -0.073074   0.050167  -1.457    0.145
## Lag2        -0.042301   0.050086  -0.845    0.398
## Lag3         0.011085   0.049939   0.222    0.824
## Lag4         0.009359   0.049974   0.187    0.851
## Lag5         0.010313   0.049511   0.208    0.835
## Volume       0.135441   0.158360   0.855    0.392
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 1731.2  on 1249  degrees of freedom
## Residual deviance: 1727.6  on 1243  degrees of freedom
## AIC: 1741.6
## 
## Number of Fisher Scoring iterations: 3

El valor p más pequeño aquí está asociado con Lag 1. El coeficiente negativo para este predictor sugiere que si el mercado tuvo un rendimiento positivo ayer, entonces es menos probable que suba hoy. Sin embargo, a un valor de 0.15, el valor p todavía es relativamente grande, por lo que no hay evidencia clara de una asociación real entre Lag 1 y Direction.

La función de predict() se puede usar para predecir la probabilidad de que el mercado suba, dados los valores de los predictores. La opción type=“response” le dice a R las probabilidades de salida. Si no se proporciona un conjunto de datos a la función predict(), entonces las probabilidades se calculan para los datos de entrenamiento que se utilizaron para ajustarse al modelo de regresión logística. Aquí hemos impreso solo las primeras diez probabilidades. Sabemos que estos valores corresponden a la probabilidad de que el mercado suba, en lugar de bajar, porque la función contrasts() indica que R ha creado una variable dummy con un 1 para Up.

glm.probs<-predict(glm.fit,type="response")
glm.probs[1:10]

##         1         2         3         4         5         6         7         8 
## 0.5070841 0.4814679 0.4811388 0.5152224 0.5107812 0.5069565 0.4926509 0.5092292 
##         9        10 
## 0.5176135 0.4888378

contrasts(Direction)

##      Up
## Down  0
## Up    1

Para hacer una predicción sobre si el mercado subirá o bajará en un día en particular, debemos convertir estas probabilidades pronosticadas en etiquetas de clase, Up o Down. Los siguientes dos comandos crean un vector de predicciones de clase en función de si la probabilidad pronosticada de un aumento del mercado es mayor o menor que 0.5

glm.pred<-rep("Down",1250)
glm.pred[glm.probs>.5]="Up"

El primer comando crea un vector de 1.250 elementos Down. La segunda línea se transforma en UP todos los elementos para los cuales la probabilidad pronosticada de un aumento de mercado excede 0.5. Dadas estas predicciones, la función table () se puede usar para producir una matriz de confusión con el fin de determinar cuántas observaciones se clasificaron correctamente o incorrectamente.

table(glm.pred,Direction)

##         Direction
## glm.pred Down  Up
##     Down  145 141
##     Up    457 507

(507+145)/1250

## [1] 0.5216

mean(glm.pred==Direction)

## [1] 0.5216

Los elementos diagonales de la matriz de confusión indican predicciones correctas, mientras que las diagonales no representan predicciones incorrectas. Por lo tanto, nuestro modelo predijo correctamente que el mercado subiría en 507 días y que bajaría en 145 días, para un total de $507 + 145 = 652$ predicciones correctas. La función mean () se puede usar para calcular la fracción de días para la cual la predicción fue correcta. En este caso, la regresión logística predijo correctamente el movimiento del mercado el 52,2% de las veces. A primera vista, parece que el modelo de regresión logística funciona un poco mejor que las conjeturas aleatorias. Sin embargo, este resultado es engañoso porque capacitamos y probamos el modelo en el mismo conjunto de 1,250 observaciones. En otras palabras, $100 − 52.2 = 47.8%$ es la tasa de terror del conjunto de entrenamiento. Como hemos visto anteriormente, la tasa de error de entrenamiento a menudo es demasiado optimista: tiende a subestimar la tasa de error de la prueba. Para evaluar mejor la precisión del modelo de regresión logística en esta configuración, podemos ajustar la parte de modelado de los datos y luego examinar qué tan bien predice los datos desactualizados. Esto generará una tasa de error más realista, en el sentido de que, en la práctica, nos interesará el rendimiento de nuestro modelo, no en los datos que solíamos ajustar al modelo, sino en los días en el futuro para los que se desconocen los movimientos del mercado. Para implementar esta estrategia, primero crearemos un vector correspondiente a las observaciones de 2001 a 2004. Luego, utilizaremos este vector para crear un conjunto de datos de observaciones de 2005.

train<-(Year<2005)#Conjunto de prueba 
Smarket.2005<-Smarket[!train,] #observaciones de 2005
dim(Smarket.2005) #252 observaciones

## [1] 252   9

Direction.2005<-Direction[!train]

El objeto traines un vector de 1.250 elementos, que corresponde a las observaciones en nuestro conjunto de datos. Los elementos del vector que corresponden a las observaciones que ocurrieron antes de 2005 se establecen en TRUE, mientras que los que corresponden a observaciones en 2005 se establecen en FALSE. El objeto train es un vector booleano, ya que sus elementos son TRUE y FAlSE. El símbolo ! Puede usarse para revertir todos los elementos de un Vector booleano Es decir, !Train es un vector similar al de entrenamiento, excepto que los elementos que son TRUE se intercambiaron a FALSE en !Train, y los elementos que son FALSE se intercambiaron a TRUE en !Train. Por lo tanto, Smarket[! Train,] produce una submatriz de los datos del mercado de valores que contiene solo las observaciones para las que se entrena FALSE, es decir, las observaciones con fechas en 2005. El resultado anterior indica que hay 252 observaciones de este tipo. Ahora ajustamos un modelo de regresión logística utilizando solo el subconjunto de las observaciones que corresponden a fechas anteriores a 2005, utilizando el argumento secundario. Luego obtenemos probabilidades pronosticadas de que el mercado de valores suba cada día de nuestro conjunto de pruebas, es decir, para los días de 2005.

glm.fit<-glm(Direction~Lag1+Lag2+Lag3+Lag4+Lag5+Volume,data=Smarket,family=binomial,subset=train)
glm.probs<-predict(glm.fit,Smarket.2005,type="response")

Ten en cuenta que hemos entrenado y probado nuestro modelo en dos conjuntos de datos completamente separados: el entrenamiento se realizó usando solo las fechas anteriores a 2005, y las pruebas se realizaron usando solo las fechas en 2005. Finalmente, calculamos las predicciones para 2005 y compárelos con los movimientos reales del mercado durante ese período de tiempo

glm.pred<-rep("Down",252)
glm.pred[glm.probs>.5]="Up"
table(glm.pred,Direction.2005)

##         Direction.2005
## glm.pred Down Up
##     Down   77 97
##     Up     34 44

mean(glm.pred==Direction.2005)

## [1] 0.4801587

mean(glm.pred!=Direction.2005)

## [1] 0.5198413

La notación != significa: No es igual a, por lo que el último comando calcula la tasa de error del conjunto de prueba. Los resultados son bastante decepcionantes: la tasa de error de la prueba es del 52%, ¡lo cual es peor que adivinar al azar! Por supuesto, este resultado no es tan sorprendente, dado que generalmente no se esperaría poder utilizar los retornos de días anteriores para predecir el rendimiento futuro del mercado. Recordamos que el modelo de regresión logística tenía valores p muy decepcionantes asociados con todos los predictores, y que el valor p más pequeño, aunque no muy pequeño, correspondía a Lag1. Quizás al eliminar las variables que parecen no ser útiles para predecir la Dirección, podemos obtener un modelo más efectivo. Después de todo, el uso de predictores que no tienen relación con la respuesta tiende a causar un deterioro en la tasa de prueba (ya que dichos predictores causan un aumento de la varianza sin una disminución correspondiente en el sesgo), por lo que eliminar dichos predictores puede producir una mejora. A continuación, hemos reajustado la regresión logística usando solo Lag1 y Lag2, que parecía tener el poder predictivo más alto en el modelo de regresión logística original.

glm.fit<-glm(Direction~Lag1+Lag2,data=Smarket,family=binomial,subset=train)
glm.probs<-predict(glm.fit,Smarket.2005,type="response")
glm.pred<-rep("Down",252)
glm.pred[glm.probs>.5]="Up"
table(glm.pred,Direction.2005)

##         Direction.2005
## glm.pred Down  Up
##     Down   35  35
##     Up     76 106

mean(glm.pred==Direction.2005)

## [1] 0.5595238

(35+106)/(252)

## [1] 0.5595238

summary (glm.fit)

## 
## Call:
## glm(formula = Direction ~ Lag1 + Lag2, family = binomial, data = Smarket, 
##     subset = train)
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)
## (Intercept)  0.03222    0.06338   0.508    0.611
## Lag1        -0.05562    0.05171  -1.076    0.282
## Lag2        -0.04449    0.05166  -0.861    0.389
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 1383.3  on 997  degrees of freedom
## Residual deviance: 1381.4  on 995  degrees of freedom
## AIC: 1387.4
## 
## Number of Fisher Scoring iterations: 3

Ahora los resultados parecen ser un poco mejores: el 56% de los movimientos diarios se han predicho correctamente. Vale la pena señalar que, en este caso, ¡una estrategia mucho más simple de predecir que el mercado aumentará todos los días también será correcta el 56% del tiempo! Por lo tanto, en términos de tasa de error general, el método de regresión logística no es mejor que el enfoque ingenuo. Sin embargo, la matriz de confusión muestra que en los días en que la regresión logística predice un aumento en el mercado, tiene una tasa de precisión del 58%. Esto sugiere una estrategia comercial posible de comprar en los días en que el modelo predice un mercado en aumento, y evitar los días de comercio en los que se predice una disminución.

3. Linear Discriminant Analysis

Ahora realizaremos LDA en los datos del Smarket. En R, ajustamos una función lda (), que es parte de la librería MASS. Observe que la sintaxis de lda () es idéntica a la de lm () y a la de glm (), excepto por la ausencia de la opción de familia. Ajustamos el modelo utilizando solo las observaciones anteriores a 2005

library(MASS)
lda.fit<-lda(Direction~Lag1+Lag2,data=Smarket,subset=train)
lda.fit

## Call:
## lda(Direction ~ Lag1 + Lag2, data = Smarket, subset = train)
## 
## Prior probabilities of groups:
##     Down       Up 
## 0.491984 0.508016 
## 
## Group means:
##             Lag1        Lag2
## Down  0.04279022  0.03389409
## Up   -0.03954635 -0.03132544
## 
## Coefficients of linear discriminants:
##             LD1
## Lag1 -0.6420190
## Lag2 -0.5135293

plot(lda.fit)

La salida LDA indica que $ˆπ1 = 0.492$ y $ˆπ2 = 0.508$; en otras palabras, el 49.2% de las observaciones de entrenamiento corresponden a días durante los cuales el mercado cayó. También proporciona las medias grupales; estos son el promedio de cada predictor dentro de cada clase, y son utilizados por LDA como estimaciones de $μk$. Esto sugiere que existe una tendencia a que los retornos de los 2 días anteriores sean negativos en los días en que el mercado aumenta, y una tendencia a que los retornos de los días anteriores sean positivos en los días en que el mercado declina. Los coeficientes de salida discriminantes lineales proporcionan la combinación lineal de Lag1 y Lag2 que se utilizan para formar la regla de decisión LDA. En otras palabras, estos son los multiplicadores de los elementos de $X = xin$ (4.19). Si $ − 0.642 × Lag1−0.514 × Lag2$ es grande, entonces el clasificador LDA predecirá un aumento del mercado, y si es pequeño, entonces el clasificador LDA predecirá una disminución del mercado. La función plot () produce gráficos de los discriminantes lineales, obtenidos calculando $− 0.642 × *Lag1* −0.514 × *Lag2*$ antes de las observaciones de entrenamiento. La función Predict () devuelve una lista con tres elementos. El primer elemento, la clase, contiene las predicciones de LDA sobre el movimiento del mercado. El segundo elemento, posterior, es una matriz cuya columna k contiene la probabilidad posterior de que la observación correspondiente pertenezca a la clase k, calculada a partir de (4.10). Finalmente, x contiene los discriminantes lineales, descritos anteriormente.

lda.pred<-predict(lda.fit, Smarket.2005)
names(lda.pred)

## [1] "class"     "posterior" "x"

Como observamos en la Sección 4.5, las predicciones de LDA y de regresión logística son casi idénticas.

lda.class<-lda.pred$class
table(lda.class,Direction.2005)

##          Direction.2005
## lda.class Down  Up
##      Down   35  35
##      Up     76 106

mean(lda.class==Direction.2005)

## [1] 0.5595238

(35+106)/252

## [1] 0.5595238

La aplicación de un umbral del 50% a las probabilidades posteriores nos permite recrear las predicciones contenidas en lda.pred$class.

sum(lda.pred$posterior[,1]>=.5)

## [1] 70

sum(lda.pred$posterior[,1]<.5)

## [1] 182

Observe que la salida de probabilidad posterior del modelo corresponde a la probabilidad de que el mercado disminuya:

lda.pred$posterior[1:20,1]

##       999      1000      1001      1002      1003      1004      1005      1006 
## 0.4901792 0.4792185 0.4668185 0.4740011 0.4927877 0.4938562 0.4951016 0.4872861 
##      1007      1008      1009      1010      1011      1012      1013      1014 
## 0.4907013 0.4844026 0.4906963 0.5119988 0.4895152 0.4706761 0.4744593 0.4799583 
##      1015      1016      1017      1018 
## 0.4935775 0.5030894 0.4978806 0.4886331

lda.class[1:20]

##  [1] Up   Up   Up   Up   Up   Up   Up   Up   Up   Up   Up   Down Up   Up   Up  
## [16] Up   Up   Down Up   Up  
## Levels: Down Up

Si quisiéramos usar un umbral de probabilidad posterior diferente del 50% para hacer predicciones, entonces podríamos hacerlo fácilmente. Por ejemplo, supongamos que deseamos predecir una disminución del mercado solo si estamos muy seguros de que el mercado disminuirá ese día, por ejemplo, si la probabilidad posterior es al menos del 90%.

sum(lda.pred$posterior[,1]>.9)

## [1] 0

Ningún día en 2005 alcanza ese umbral! De hecho, la mayor probabilidad de disminución posterior en todo 2005 fue de 52.02%

4. Quadratic Discriminant Analysis

Ahora ajustaremos un modelo QDA a los datos de Smarket. QDA se implementa al utilizar la función qda (), que también forma parte de la librería MASS. La sintaxis de qda () es idéntica a la de lda ().

qda.fit<-qda(Direction~Lag1+Lag2,data=Smarket,subset=train)
qda.fit

## Call:
## qda(Direction ~ Lag1 + Lag2, data = Smarket, subset = train)
## 
## Prior probabilities of groups:
##     Down       Up 
## 0.491984 0.508016 
## 
## Group means:
##             Lag1        Lag2
## Down  0.04279022  0.03389409
## Up   -0.03954635 -0.03132544

La salida contiene las medias del grupo. Pero no contiene los coeficientes de los discriminantes lineales, porque el clasificador QDA implica una función de predictores cuadráticos, más que lineal. La función predict () funciona exactamente de la misma manera que para LDA.

qda.class<-predict(qda.fit,Smarket.2005)$class
table(qda.class,Direction.2005)

##          Direction.2005
## qda.class Down  Up
##      Down   30  20
##      Up     81 121

mean(qda.class==Direction.2005)

## [1] 0.5992063

Curiosamente, las predicciones de QDA son precisas casi el 60% del tiempo, a pesar de que los datos de 2005 no se utilizaron para ajustarse al modelo. Este nivel de precisión es bastante impresionante para los datos del mercado de valores, que se sabe que es bastante difícil de modelar con precisión. Esto sugiere que la forma cuadrática asumida por QDA puede capturar la relación verdadera con mayor precisión que las formas lineales asumidas por LDA y la regresión logística. Sin embargo, recomendamos evaluar el rendimiento de este método en un conjunto de prueba más grande antes de apostar a que este enfoque superará constantemente al mercado.

5. K-Nearest Neighbors

Ahora realizaremos KNN usando la función knn (), que es parte de la libreria class. Esta función funciona de manera bastante diferente de las otras funciones de ajuste del modelo que hemos encontrado hasta ahora. En lugar de un enfoque de dos pasos en el que primero ajustamos el modelo y luego usamos el modelo para hacer predicciones, knn () forma predicciones usando un solo comando. La función requiere cuatro entradas.

Una matriz que contiene los predictores asociados con los datos de entrenamiento, etiquetados como train.X a continuación.
Una matriz que contiene los predictores asociados con los datos para los que deseamos hacer predicciones, etiquetada como test.X a continuación.
Un vector que contiene las etiquetas de clase para las observaciones de entrenamiento, etiquetado train.Direction
Un valor para K, el número de vecinos más cercanos que utilizará el clasificador.

Usamos la función cbind (), abreviatura de enlace de columna, para unir las variables Lag1 y Lag2 juntas en dos matrices, una para el conjunto de entrenamiento y la otra para el conjunto de prueba.

library(class)
train.X<-cbind(Lag1,Lag2)[train,]
test.X<-cbind(Lag1,Lag2)[!train,]
train.Direction<-Direction[train]

Ahora la función knn () se puede usar para predecir el movimiento del mercado para las fechas en 2005. Establecemos una semilla aleatoria antes de aplicar knn () porque si varias observaciones están vinculadas como vecinas más cercanas, entonces Rompera el empate al azar. Por lo tanto, una semilla debe asediar para garantizar la reproducibilidad de los resultados.

set.seed(1)
knn.pred<-knn(train.X,test.X,train.Direction,k=1)
table(knn.pred,Direction.2005)

##         Direction.2005
## knn.pred Down Up
##     Down   43 58
##     Up     68 83

(83+43)/252

## [1] 0.5

Los resultados con K = 1 no son muy buenos, ya que solo el 50% de las observaciones se predicen correctamente. Por supuesto, puede ser que K = 1 resulte en un ajuste demasiado flexible a los datos. A continuación, repetimos el análisis con K = 3.

knn.pred<-knn(train.X,test.X,train.Direction,k=3)
table(knn.pred,Direction.2005)

##         Direction.2005
## knn.pred Down Up
##     Down   48 54
##     Up     63 87

mean(knn.pred==Direction.2005)

## [1] 0.5357143

Los resultados han mejorado ligeramente. Pero el aumento de K más lejos resulta no proporcionar más mejoras. Parece que para estos datos, QDA proporciona los mejores resultados de los métodos que hemos examinado hasta ahora.

6. Resuelve

Utilizando el conjunto de datos de Boston, ajuste los modelos de clasificación para predecir si un suburbio determinado tiene una tasa de criminalidad superior o inferior a la mediana. Explore el modelo de regresión logística utilizando varios subconjuntos de predictores. Describe tus hallazgos.

# Minimal Boston Logistic Classification
# Requiere solo MASS (glm es de stats base)
library(MASS)
set.seed(123)

# 1) Variable objetivo binaria (HighCrime) y split train/test
bdat <- Boston
bdat$HighCrime <- factor(ifelse(bdat$crim > median(bdat$crim), "High", "Low"),
                         levels = c("Low", "High"))
idx <- sample(nrow(bdat), round(0.7 * nrow(bdat)))
train <- bdat[idx, ]
test  <- bdat[-idx, ]

En esta prte caragamos el paquete de Mass donde vienen los datos de Boston y creamos la variable objetivo binaria HighCrime que nos indica si la tasa de criminalidad es alta o baja en comparación con la mediana. Luego, dividimos los datos en conjuntos de entrenamiento y prueba, utilizando el 70% de los datos para el entrenamiento y el 30% restante para la prueba.

# 2) Fórmulas: completo (sin 'crim') y subconjuntos
forms <- list(
  full    = HighCrime ~ . - crim,
  socio   = HighCrime ~ lstat + medv + ptratio + black,
  estruct = HighCrime ~ nox + rm + age + dis + chas,
  acceso  = HighCrime ~ rad + tax + indus + zn
)

Aquí definimos varias fórmulas para diferentes modelos de regresión logística. El modelo completo utiliza todas las variables predictoras excepto ‘crim’, mientras que los otros modelos utilizan subconjuntos específicos de variables relacionadas con factores socioeconómicos, estructurales y de acceso.

# 3) Función corta de ajuste y evaluación (umbral 0.5)
fit_eval <- function(fm) {
  m <- glm(fm, data = train, family = binomial)
  p <- predict(m, newdata = test, type = "response")
  pred <- factor(ifelse(p >= 0.5, "High", "Low"), levels = levels(test$HighCrime))
  acc <- mean(pred == test$HighCrime)
  list(acc = acc, cm = table(Pred = pred, Obs = test$HighCrime))
}

Aquí definimos una función llamada fit_eval que ajusta un modelo de regresión logística utilizando una fórmula dada, realiza predicciones en el conjunto de prueba y evalúa la precisión del modelo. La función devuelve la precisión y la matriz de confusión.

# 4) Resultados por modelo
res <- lapply(forms, fit_eval)

# Accuracy por modelo
sapply(res, `[[`, "acc")

##      full     socio   estruct    acceso 
## 0.8881579 0.7565789 0.8750000 0.8486842

# Matrices de confusión
res$full$cm

##       Obs
## Pred   Low High
##   Low   64    6
##   High  11   71

res$socio$cm

##       Obs
## Pred   Low High
##   Low   57   19
##   High  18   58

res$estruct$cm

##       Obs
## Pred   Low High
##   Low   62    6
##   High  13   71

res$acceso$cm

##       Obs
## Pred   Low High
##   Low   67   15
##   High   8   62

summary(glm(HighCrime ~ . - crim, data=train, family=binomial))

## 
## Call:
## glm(formula = HighCrime ~ . - crim, family = binomial, data = train)
## 
## Coefficients:
##               Estimate Std. Error z value Pr(>|z|)    
## (Intercept) -28.064366   8.944009  -3.138  0.00170 ** 
## zn           -0.093211   0.042125  -2.213  0.02692 *  
## indus        -0.044137   0.058791  -0.751  0.45281    
## chas          0.247705   0.823661   0.301  0.76362    
## nox          42.153836   9.034180   4.666 3.07e-06 ***
## rm            0.033207   0.882360   0.038  0.96998    
## age           0.028651   0.014189   2.019  0.04346 *  
## dis           0.657709   0.258111   2.548  0.01083 *  
## rad           0.532624   0.183178   2.908  0.00364 ** 
## tax          -0.007883   0.003367  -2.341  0.01923 *  
## ptratio       0.448510   0.161183   2.783  0.00539 ** 
## black        -0.034710   0.014439  -2.404  0.01622 *  
## lstat         0.103425   0.061650   1.678  0.09342 .  
## medv          0.220912   0.090027   2.454  0.01413 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 490.74  on 353  degrees of freedom
## Residual deviance: 145.64  on 340  degrees of freedom
## AIC: 173.64
## 
## Number of Fisher Scoring iterations: 9

El modelo Full tiene mayor exactitud, pero el modelo socioeconómico no está muy lejos. Los modelos estructurales y de acceso tienen un rendimiento significativamente peor. Esto sugiere que los factores socioeconómicos son los más importantes para predecir la tasa de criminalidad en los suburbios de Boston. El análisis de los coeficientes del modelo completo revela que variables como lstat (porcentaje de población de estatus bajo), medv (valor medio de las viviendas) y ptratio (relación alumno-maestro) son estadísticamente significativas, lo que respalda la importancia de los factores socioeconómicos en la predicción de la criminalidad.