Ćwiczenie 2.2 elementów uczenia statystycznego

Podręcznik najpierw generuje dane 2-klasowe poprzez:

co daje:

a następnie pyta:

Próbuję rozwiązać ten problem, najpierw modelując to za pomocą tego modelu graficznego:

gdzie $c$ jest etykieta, $h\,(1\le h \le 10)$ jest indeksem wybranej średniej $m_h^c$, i $x$jest punktem danych. To da

$$\begin{align*} \Pr(x\mid m_h^c) =& \mathcal{N}(m_h^c,\mathbf{I}/5)\\ \Pr(m_h^c\mid h,c=\mathrm{blue}) =& \mathcal{N}((1,0)^T,\mathbf{I})\\ \Pr(m_h^c\mid h,c=\mathrm{orange}) =& \mathcal{N}((0,1)^T,\mathbf{I})\\ \Pr(h) =& \frac{1}{10}\\ \Pr(c) =& \frac{1}{2} \end{align*}$$

Z drugiej strony granicą jest . Dzięki rządowi bayesowskiemu mamy$\{x:\Pr(c=\mathrm{blue}\mid x)=\Pr(c=\mathrm{orange}\mid x)\}$

$$\begin{align*} \Pr(c\mid x) =& \frac{\Pr(x\mid c)\Pr(c)}{\sum_c\Pr(x\mid c)\Pr(c)}\\ \Pr(x\mid c) =& \sum_h\int_{m_h^c}\Pr(h)\Pr(m_h^c\mid h,c)\Pr(x\mid m_h^c) \end{align*}$$

Ale później odkryłem, że ustawienie problemu jest symetryczne, więc może to dać jako granicę. Jeśli problemem jest ustalenie granicy, kiedy są uwarunkowane, równanie obejmie parametrów, które moim zdaniem raczej nie będą celem ćwiczenia.$x=y$$m_h^c$$40$

Czy coś rozumiem? Dziękuję Ci.

Nie sądzę, że powinieneś znaleźć wyrażenie analityczne dla granicy decyzji Bayesa, dla danej realizacji . Podobnie wątpię, abyś miał granicę nad rozkładem , ponieważ to tylko symetrii, jak zauważyłeś.$m_k$$m_k$$x=y$

Myślę, że musisz pokazać program, który może obliczyć granicę decyzji dla danej realizacji . Można to zrobić poprzez ustawienie w dół siatkę i wartości, obliczanie gęstości klasy-warunkowy, a znalezienie punktów, gdzie są one równe.$m_k$$x$$y$

Ten kod to dźgnięcie. IIRC faktycznie ma kod do obliczania granicy decyzji w Modern Applied Statistics z S , ale nie mam tego teraz przydatnego.

# for dmvnorm/rmvnorm: multivariate normal distribution
library(mvtnorm)

# class-conditional density given mixture centers
f <- function(x, m)
{
    out <- numeric(nrow(x))
    for(i in seq_len(nrow(m)))
        out <- out + dmvnorm(x, m[i, ], diag(0.2, 2))
    out
}

# generate the class mixture centers
m1 <- rmvnorm(10, c(1,0), diag(2))
m2 <- rmvnorm(10, c(0,1), diag(2))
# and plot them
plot(m1, xlim=c(-2, 3), ylim=c(-2, 3), col="blue")
points(m2, col="red")

# display contours of the class-conditional densities
dens <- local({
    x <- y <- seq(-3, 4, len=701)
    f1 <- outer(x, y, function(x, y) f(cbind(x, y), m1))
    f2 <- outer(x, y, function(x, y) f(cbind(x, y), m2))
    list(x=x, y=y, f1=f1, f2=f2)
})

contour(dens$x, dens$y, dens$f1, col="lightblue", lty=2, levels=seq(.3, 3, len=10),
        labels="", add=TRUE)

contour(dens$x, dens$y, dens$f2, col="pink", lty=2, levels=seq(.3, 3, len=10),
        labels="", add=TRUE)

# find which points are on the Bayes decision boundary
eq <- local({
    f1 <- dens$f1
    f2 <- dens$f2
    pts <- seq(-3, 4, len=701)
    eq <- which(abs((dens$f1 - dens$f2)/(dens$f1 + dens$f2)) < 5e-3, arr.ind=TRUE)
    eq[,1] <- pts[eq[,1]]
    eq[,2] <- pts[eq[,2]]
    eq
})
points(eq, pch=16, cex=0.5, col="grey")

---

Wynik:

W rzeczywistości książka wymaga analitycznego rozwiązania tego problemu. I tak, musisz uwarunkować granicę, ale nie na 40: nigdy nie poznasz ich dokładnie. Zamiast tego musisz uzależnić 200 punktów danych, które zobaczysz. Będziesz potrzebował 200 parametrów, ale ze względu na użycie sumowania odpowiedź nie wydaje się zbyt skomplikowana.

Nigdy nie byłbym w stanie wyprowadzić tej formuły, więc biorę jedynie za to, że uświadomiłem sobie, że rozwiązanie analityczne nie musi być brzydkie, a następnie szukałem go w Google. Na szczęście, to świadczone przez autorów, niektórzy mili ludzie, strony 6-7 .

Chciałbym wcześniej natknąć się na powyższy kod; właśnie zakończyłem tworzenie alternatywnego kodu poniżej ... za to, co jest warte

set.seed(1)
library(MASS)

#create original 10 center points/means for each class 
I.mat=diag(2)
mu1=c(1,0);mu2=c(0,1)
mv.dist1=mvrnorm(n = 10, mu1, I.mat)
mv.dist2=mvrnorm(n = 10, mu2, I.mat)

values1=NULL;values2=NULL

#create 100 observations for each class, after random sampling of a center point, based on an assumed bivariate probability distribution around each center point  
for(i in 1:10){
  mv.values1=mv.dist1[sample(nrow(mv.dist1),size=1,replace=TRUE),]
  sub.mv.dist1=mvrnorm(n = 10, mv.values1, I.mat/5)
  values1=rbind(sub.mv.dist1,values1)
}
values1

#similar as per above, for second class
for(i in 1:10){
  mv.values2=mv.dist2[sample(nrow(mv.dist2),size=1,replace=TRUE),]
  sub.mv.dist2=mvrnorm(n = 10, mv.values2, I.mat/5)
  values2=rbind(sub.mv.dist2,values2)
}
values2

#did not find probability function in MASS, so used mnormt
library(mnormt)

#create grid of points
grid.vector1=seq(-2,2,0.1)
grid.vector2=seq(-2,2,0.1)
length(grid.vector1)*length(grid.vector2)
grid=expand.grid(grid.vector1,grid.vector2)



#calculate density for each point on grid for each of the 100 multivariates distributions
prob.1=matrix(0:0,nrow=1681,ncol=10) #initialize grid
for (i in 1:1681){
  for (j in 1:10){
    prob.1[i,j]=dmnorm(grid[i,], mv.dist1[j,], I.mat/5)  
  }
}
prob.1
prob1.max=apply(prob.1,1,max)

#second class - as per above
prob.2=matrix(0:0,nrow=1681,ncol=10) #initialize grid
for (i in 1:1681){
  for (j in 1:10){
    prob.2[i,j]=dmnorm(grid[i,], mv.dist2[j,], I.mat/5)  
  }
}
prob.2
prob2.max=apply(prob.2,1,max)

#bind
prob.total=cbind(prob1.max,prob2.max)
class=rep(1,1681)
class[prob1.max<prob2.max]=2
cbind(prob.total,class)

#plot points
plot(grid[,1], grid[,2],pch=".", cex=3,col=ifelse(class==1, "coral", "cornflowerblue"))

points(values1,col="coral")
points(values2,col="cornflowerblue")

#check - original centers
# points(mv.dist1,col="coral")
# points(mv.dist2,col="cornflowerblue")