Para este ejercisio se utilzó la base de datos seguros.csv. El objetivo consiste en construir una red bayesiana para la cual se deben considerar algunas reglas propuestas por un experto en seguros y reglas propuestas por el usuario en caso de ser necesario.
El experto en seguros ha establecido para la construcción e interpretación de la red, las siguientes reglas:
La base de datos seguros.csv esta compuesta de la siguiente forma
seguros <- read.csv("seguros.csv")
attach(seguros)
head(seguros)
## Age SocioEcon RiskAversion VehicleYear Accident MakeModel
## 14494 Adult Middle Normal Current Severe Economy
## 18950 Senior Prole Normal Older None SportsCar
## 13705 Adult UpperMiddle Cautious Current None Luxury
## 17328 Senior Prole Cautious Older None Economy
## 19036 Adult UpperMiddle Normal Current None SportsCar
## 16971 Senior Prole Normal Older None Economy
## DrivQuality Antilock CarValue MedCost Airbag DrivHist
## 14494 Poor False TenThou Thousand True Many
## 18950 Normal False TenThou Thousand False Zero
## 13705 Excellent True FiftyThou Thousand True Zero
## 17328 Normal False FiveThou Thousand False Zero
## 19036 Normal True TenThou Thousand True Zero
## 16971 Normal False FiveThou Thousand False One
Para la construcción de la red de tal manera que considere las restricciones impuestas por el experto en seguros, es necesario crear una blacklist que prohiba la existencia de arcos entre ciertos nodos y usar la función hill climbing para encontrar una posible estructura de la red a partir de los datos.
library("bnlearn")
black <- data.frame(from = c("Airbag", "Airbag", "Airbag", "Accident", "Accident",
"Age", "SocioEcon", "RiskAversion", "Accident", "RiskAversion", "DrivQuality"),
to = c("DrivQuality", "DrivHist", "Accident", "MakeModel", "VehicleYear",
"MedCost", "MedCost", "SocioEcon", "SocioEcon", "Age", "Age"))
net.seguros <- hc(seguros, score = "aic", blacklist = black)
net.seguros
##
## Bayesian network learned via Score-based methods
##
## model:
## [Age][SocioEcon|Age][VehicleYear|SocioEcon][MakeModel|SocioEcon]
## [Antilock|VehicleYear:MakeModel][CarValue|VehicleYear:MakeModel]
## [Airbag|VehicleYear:MakeModel][Accident|Age:Antilock]
## [DrivQuality|Age:Accident][MedCost|Accident:MakeModel]
## [RiskAversion|Age:SocioEcon:DrivQuality]
## [DrivHist|RiskAversion:DrivQuality]
## nodes: 12
## arcs: 20
## undirected arcs: 0
## directed arcs: 20
## average markov blanket size: 4.00
## average neighbourhood size: 3.33
## average branching factor: 1.67
##
## learning algorithm: Hill-Climbing
## score: AIC (disc.)
## penalization coefficient: 1
## tests used in the learning procedure: 314
## optimized: TRUE
La estructura que se obtuvo es la siguiente:
graphviz.plot(net.seguros)
La misma estructura representada de otra forma esta dada por:
plot(net.seguros)
Una vez incorporada la información del experto en la estimación de la red, vale la pena incorporar dos restricciones adiconales, debido a que es razonable pesnar que el historial de manejo DrivHist de una persona debe depender de la edad de la persona y que el hecho de que un carro tenga Antilock no determina si puede o no tener un accidente, por esta razón al incluir estas nuevas restricciones la red bayesiana que se obtiene es la siguiente:
black <- data.frame(from = c("Airbag", "Airbag", "Airbag", "Accident", "Accident",
"Age", "SocioEcon", "RiskAversion", "Accident", "RiskAversion", "DrivQuality",
"Antilock"), to = c("DrivQuality", "DrivHist", "Accident", "MakeModel",
"VehicleYear", "MedCost", "MedCost", "SocioEcon", "SocioEcon", "Age", "Age",
"Accident"))
white <- data.frame(from = c("Age"), to = c("DrivHist"))
net.seguros <- hc(seguros, score = "aic", blacklist = black, whitelist = white)
net.seguros
##
## Bayesian network learned via Score-based methods
##
## model:
## [Age][SocioEcon|Age][RiskAversion|Age:SocioEcon][VehicleYear|SocioEcon]
## [MakeModel|SocioEcon][Antilock|VehicleYear:MakeModel]
## [CarValue|VehicleYear:MakeModel][Airbag|VehicleYear:MakeModel]
## [DrivHist|Age:RiskAversion][DrivQuality|Age:RiskAversion:DrivHist]
## [Accident|DrivQuality:CarValue][MedCost|Accident:MakeModel]
## nodes: 12
## arcs: 20
## undirected arcs: 0
## directed arcs: 20
## average markov blanket size: 3.83
## average neighbourhood size: 3.33
## average branching factor: 1.67
##
## learning algorithm: Hill-Climbing
## score: AIC (disc.)
## penalization coefficient: 1
## tests used in the learning procedure: 332
## optimized: TRUE
cuya gráfica es
graphviz.plot(net.seguros)
Ya que se ha definido una posible estructura de la red que represente los datos, es factible pensar en la construcción de modelos locales que permitan tener mejores estimaciones de las probabilidades condicionales. Para ello, se propone usar un modelo logístico multinomial para describir la relación VehicleYear-> CarValue <- MakeModel pues la variable CarValue tiene más de dos niveles. Para evaluar la precisión del modelo se ha dividio la base en muestra de entrenamiento y muestra de prueba. Ajustando el modelo con los datos pertencientes a la muestra de entrenamiento
library(glmnet)
muestra.reng <- sample(1:nrow(seguros), 3500)
seguros.muestra <- seguros[muestra.reng, ]
mat.1 <- model.matrix(~-1 + VehicleYear + MakeModel + VehicleYear:MakeModel,
data = seguros.muestra[, c("VehicleYear", "MakeModel")])
mod.cv <- cv.glmnet(y = seguros.muestra$CarValue, x = mat.1, alpha = 0.6, family = "multinomial")
plot(mod.cv)
Hacemos las predicciones para los datos originales:
mat.2 <- model.matrix(~-1 + VehicleYear + MakeModel + VehicleYear:MakeModel,
data = seguros[, c("VehicleYear", "MakeModel")])
mod.cv.pred <- predict(mod.cv, s = mod.cv$lambda.min, type = "response", newx = mat.2)[,
, 1]
head(mod.cv.pred)
## FiftyThou FiveThou Million TenThou TwentyThou
## 14494 7.762e-04 8.643e-02 1.972e-04 0.8033296 0.1092671
## 18950 2.889e-02 2.743e-01 1.086e-02 0.4006369 0.2853148
## 13705 9.994e-01 9.273e-05 4.100e-07 0.0003306 0.0002199
## 17328 3.366e-05 7.845e-01 2.138e-05 0.2063452 0.0090708
## 19036 1.262e-01 1.154e-03 1.897e-02 0.1174842 0.7362064
## 16971 3.366e-05 7.845e-01 2.138e-05 0.2063452 0.0090708
Dado que a estructura es un colisionador, las variablesVehicleYear y MakeModel son independientes a menos que se tenga información del valor de CarValue, por ejemplo, si se conoce que el valor del coche es de TenThou, las predicciones del modelo son:
library(ggplot2)
library(plyr)
seguros$pred.1 <- mod.cv.pred[, 4]
dat.res <- ddply(seguros, c("VehicleYear", "MakeModel"), summarise, prob = mean(pred.1))
ggplot(dat.res, aes(x = MakeModel, y = prob, colour = VehicleYear, group = VehicleYear)) +
geom_point() + geom_line()
lo cual nos muestra que dado que es un carro con valor de TenThou las probabilidades de que se trate de un carro actual van disminuyendo conforme aumenta la clase del auto, además de que no se observan carros de super lujo con ese precio.
De igual forma podemos realizar el mismo análsis para autos de valor TwentyThou las probabilidades de que se trate de un carro viejo aumentan conforme aumenta la categoría del automovil y dado ese valor, prácticamente solo encontraremos carros del tipo FamilySedan y SportsCar que sean actuales.
library(ggplot2)
library(plyr)
seguros$pred.1 <- mod.cv.pred[, 5]
dat.res <- ddply(seguros, c("VehicleYear", "MakeModel"), summarise, prob = mean(pred.1))
ggplot(dat.res, aes(x = MakeModel, y = prob, colour = VehicleYear, group = VehicleYear)) +
geom_point() + geom_line()
