

## EM for simple regression

## data

x<-rnorm(100)
y<-x+.5*rnorm(100)

summary(lm(y~x-1))

y[91:100]<-NA

a0<-1
s0<-.5

e<-1

while (e>.00001){

sx<-sum(x^2)
sy<-sum(y[1:90]^2)+a0^2*sum(x[91:100]^2)+10*s0
sxy<-sum(x[1:90]*y[1:90])+a0*sum(x[91:100]^2)

a<-sxy/sx
s<-(sx*sy-sxy^2)/(sx^2)

e<-(a0-a)^2+(s0-s)^2

a0<-a
s0<-s
}