Comment tracer la frontière de décision d'un classificateur k-plus proche voisin à partir des éléments d'apprentissage statistique?

31

Je veux générer l'intrigue décrite dans le livre ElemStatLearn "The Elements of Statistical Learning: Data Mining, Inference, and Prediction. Second Edition" de Trevor Hastie & Robert Tibshirani & Jerome Friedman. L'intrigue est:

entrez la description de l'image ici

Je me demande comment je peux produire ce graphique exact dans R, notez en particulier les graphiques de la grille et le calcul pour montrer la frontière.

littleEinstein
la source
1
@StasK: oui, ça l'est. Comment générer l'intrigue? Pourriez-vous s'il vous plaît aider? Merci beaucoup!
littleEinstein

Réponses:

35

Pour reproduire cette figure, vous devez avoir installé le package ElemStatLearn sur votre système. L'ensemble de données artificielles a été généré avec mixture.example()comme indiqué par @StasK.

library(ElemStatLearn)
require(class)
x <- mixture.example$x
g <- mixture.example$y
xnew <- mixture.example$xnew
mod15 <- knn(x, xnew, g, k=15, prob=TRUE)
prob <- attr(mod15, "prob")
prob <- ifelse(mod15=="1", prob, 1-prob)
px1 <- mixture.example$px1
px2 <- mixture.example$px2
prob15 <- matrix(prob, length(px1), length(px2))
par(mar=rep(2,4))
contour(px1, px2, prob15, levels=0.5, labels="", xlab="", ylab="", main=
        "15-nearest neighbour", axes=FALSE)
points(x, col=ifelse(g==1, "coral", "cornflowerblue"))
gd <- expand.grid(x=px1, y=px2)
points(gd, pch=".", cex=1.2, col=ifelse(prob15>0.5, "coral", "cornflowerblue"))
box()

Toutes les commandes, à l'exception des trois dernières, proviennent de l'aide en ligne de mixture.example. Notez que nous avons utilisé le fait que expand.gridnous organiserons sa sortie en variant d' xabord, ce qui permet en outre d'indexer (par colonne) les couleurs dans la prob15matrice (de dimension 69x99), qui contient la proportion des votes pour la classe gagnante pour chaque coordonnée de réseau ( px1, px2).

entrez la description de l'image ici

chl
la source
+1. Merci! Je me demande également comment générer les données comme décrit dans le texte "exposer l'oracle". Pourriez-vous également ajouter cela au lieu d'utiliser les données du site Web?
littleEinstein
@littleEinstein Voulez-vous dire à quoi sert l'aide en ligne mixture.example? Regardez la configuration de simulation sous la ligne commençant par # Reproducing figure 2.4, page 17 of the book:dans la section exemple.
chl
pouvez-vous s'il vous plaît laissez-moi savoir le lien? Je ne le trouve pas.
littleEinstein
Désolé @littleEinstein, mais il me manque probablement quelque chose. Il s'agit simplement de taper help(mixture.example)ou example(mixture.example)à l'invite R (après avoir chargé le package requis avec library(ElemStatLearn)). Le code pour générer le jeu de données artificiel (pour ne pas générer la figure 2.4) est écrit en clair R dans la section Exemple.
chl
1
BTW, je viens de tomber sur le blog de @ Shane où il a utilisé ggplotà des fins similaires. Vérifiez ceci: ESL 2.1: régression linéaire vs KNN .
chl
7

Je suis auto-apprenant ESL et j'essaie de travailler à travers tous les exemples fournis dans le livre. Je viens de le faire et vous pouvez vérifier le code R ci-dessous:

library(MASS)
# set the seed to reproduce data generation in the future
seed <- 123456
set.seed(seed)

# generate two classes means
Sigma <- matrix(c(1,0,0,1),nrow = 2, ncol = 2)
means_1 <- mvrnorm(n = 10, mu = c(1,0), Sigma)
means_2 <- mvrnorm(n = 10, mu = c(0,1), Sigma)

# pick an m_k at random with probability 1/10
# function to generate observations
genObs <- function(classMean, classSigma, size, ...)
{
  # check input
  if(!is.matrix(classMean)) stop("classMean should be a matrix")
  nc <- ncol(classMean)
  nr <- nrow(classMean)
  if(nc != 2) stop("classMean should be a matrix with 2 columns")
  if(ncol(classSigma) != 2) stop("the dimension of classSigma is wrong")

  # mean for each obs
    # pick an m_k at random
  meanObs <- classMean[sample(1:nr, size = size, replace = TRUE),]
  obs <- t(apply(meanObs, 1, function(x) mvrnorm(n = 1, mu = x, Sigma = classSigma )) )
  colnames(obs) <- c('x1','x2')
  return(obs)
}


obs100_1 <- genObs(classMean = means_1, classSigma = Sigma/5, size = 100)
obs100_2 <- genObs(classMean = means_2, classSigma = Sigma/5, size = 100)

# generate label
y <- rep(c(0,1), each = 100)

# training data matrix
trainMat <- as.data.frame(cbind(y, rbind(obs100_1, obs100_2)))

# plot them
library(lattice)
with(trainMat, xyplot(x2 ~ x1,groups = y, col=c('blue', 'orange')))

# now fit two models

# model 1: linear regression
lmfits <- lm(y ~ x1 + x2 , data = trainMat)

# get the slope and intercept for the decision boundary
intercept <- -(lmfits$coef[1] - 0.5) / lmfits$coef[3]
slope <- - lmfits$coef[2] / lmfits$coef[3]

# Figure 2.1
xyplot(x2 ~ x1, groups = y, col = c('blue', 'orange'), data = trainMat,
       panel = function(...)
       {
        panel.xyplot(...)
        panel.abline(intercept, slope)
        },
       main = 'Linear Regression of 0/1 Response')    

# model2: k nearest-neighbor methods
library(class)
# get the range of x1 and x2
rx1 <- range(trainMat$x1)
rx2 <- range(trainMat$x2)
# get lattice points in predictor space
px1 <- seq(from = rx1[1], to = rx1[2], by = 0.1 )
px2 <- seq(from = rx2[1], to = rx2[2], by = 0.1 )
xnew <- expand.grid(x1 = px1, x2 = px2)

# get the contour map
knn15 <- knn(train = trainMat[,2:3], test = xnew, cl = trainMat[,1], k = 15, prob = TRUE)
prob <- attr(knn15, "prob")
prob <- ifelse(knn15=="1", prob, 1-prob)
prob15 <- matrix(prob, nrow = length(px1), ncol = length(px2))

# Figure 2.2
par(mar = rep(2,4))
contour(px1, px2, prob15, levels=0.5, labels="", xlab="", ylab="", main=
    "15-nearest neighbour", axes=FALSE)
points(trainMat[,2:3], col=ifelse(trainMat[,1]==1, "coral", "cornflowerblue"))
points(xnew, pch=".", cex=1.2, col=ifelse(prob15>0.5, "coral", "cornflowerblue"))
box()
Daoying Lin
la source
1
Pour entrer du code ici sans faire cela, vous pouvez mettre en surbrillance le texte qui est du code puis cliquer sur le bouton "code" en haut de la page. C'est dans une rangée d'icônes / boutons. Le code ressemble à des accolades.
Peter Flom - Réintègre Monica
Re: "comment coller un bloc de code R". Vous avez accès à une petite barre de menus lors de la modification de votre message.
chl
De plus, si vous n'utilisez pas un éditeur qui peut indenter facilement des blocs de code, je pense que vous serez heureux de passer à un éditeur. Par exemple, dans Rstudio, sélectionner le code et appuyer sur tab le met en retrait, dans vim vous pouvez 5>>, etc.
Mark