# Konvoluutio 19.1.2011
########################

x<-0:10
# Poissonin jakauma Poi(0.5)
plot(x,dpois(x,0.5),type="h",xlab="x",ylab="f(x)",xlim=c(0,10),ylim=c(0,0.61))
title("Poisson(0.5)")

# Binomijakauma Bin(5,0.3)
plot(x,dbinom(x,5,0.3),type="h",xlab="x",ylab="f(v)",xlim=c(0,10),ylim=c(0,0.4))
title("binomi(5,0.3)")

# Vierekkäin Poissonin jakauma ja Binomijakauma
par(mfrow=c(1,2))
plot(x,dpois(x,0.5),type="h",xlab="x",ylab="f(x)",xlim=c(0,10),ylim=c(0,0.61))
title("Poisson(0.5)")
plot(x,dbinom(x,5,0.3),type="h",xlab="x",ylab="f(v)",xlim=c(0,10),ylim=c(0,0.4))
title("binomi(5,0.3)")

# Y on Poisson ja Binomi-satunnaismuuttujien summa
# Esim. todennäköisyys P(Y=y), kun y=3 
y<-3
dpois(0:y,0.5)
dbinom(y:0,5,0.3)

# Konvoluutio
dpois(0:y,0.5)*dbinom(y:0,5,0.3)
sum(dpois(0:y,0.5)*dbinom(y:0,5,0.3))

for (y in 0:15) print(sum(dpois(0:y,0.5)*dbinom(y:0,5,0.3)))

yr<-15
fy<-rep(0,yr)
for (y in 0:yr){
fy[y+1]<-sum(dpois(0:y,0.5)*dbinom(y:0,5,0.3))
}
# fy[y+1]=p(Y=y)

y<-0:yr
plot(y,fy,type="h",xlab="x",ylab="f(y)",xlim=c(0,10),ylim=c(0,0.4))
title("Summan jakauma")

par(mfrow=c(2,2))
x<-0:10
plot(x,dpois(x,0.5),type="h",xlab="x",ylab="f(x)",xlim=c(0,10),ylim=c(0,0.61))
title("Poisson(0.5)")
plot(x,dbinom(x,5,0.3),type="h",xlab="x",ylab="f(v)",xlim=c(0,10),ylim=c(0,0.4))
title("binomi(5,0.3)")
y<-0:yr
plot(y,fy,type="h",xlab="x",ylab="f(y)",xlim=c(0,10),ylim=c(0,0.4))
title("Summan jakauma")

#####################
konvoluutio<-function(yr){
fy<-rep(0,yr)
for (y in 0:yr){
fy[y+1]<-sum(dpois(0:y,0.5)*dbinom(y:0,5,0.3))
}
return(fy)
}

konvoluutio(10)
y<-0:10
plot(y,konvoluutio(10),type="hist")