# normal data, in 3d
library(emdbook)
curve3d( dnorm( y , mean=1 + 1/2*x , sd=2 ) , from=c(-15,-10) , to=c(15,10) , n=c(40,40) , theta=-50 , phi=40 , zlab="probability" )

# another way to view it
plot( c(0,1) ~ c(0,1) , type="l" , xlab="x" , ylab="y" )
y <- rnorm( 10000 , mean=0 , sd=2 )
hist(y,axes=FALSE,xlab="",ylab="",main="")

# parabola example
x <- rnorm( 30 , 9 , 4 )
y <- rnorm( 30 , 10 + 4*x^2 , 30 )
plot( y ~ x )

x2 <- x^2
summary( lm( y ~ x2 ) )

# exponential example
x <- rnorm( 30 , 4 , 3 )
y <- rnorm( 30 , mean=2*exp( -1/2*x ) , 0.15 )
plot( y[x>0] ~ x[x>0] )

# fit
library(bbmle)
m <- mle2( y ~ dnorm( a*exp(b*x) , s ) , start=list(a=2,b=-0.5,s=0.1) )

plot( y[x>0] ~ x[x>0] )
curve( coef(m)[1] * exp( coef(m)[2] * x ) , from=0 , to=12 , add=TRUE )

# binomial distribution examples

# utility function
prep.bar <- function( d , dmax=max(d) ) {
    x <- rep(0,dmax+1)
    for ( i in 1:(dmax+1) ) {
        x[i] <- length( d[d==(i-1)] )
    }
    x
}

y <- rbinom( 10000 , prob=0.1 , size=10 )
barplot( prep.bar(y) , names.arg=0:max(y) , xlab="y" , ylab="Frequency" )

# show likelihood imposed over barplot
ss <- 100
y <- rbinom( ss , prob=0.1 , size=10 )
barplot( prep.bar(y) , names.arg=0:max(y) , xlab="y" , ylab="Frequency" , space=0 )
p <- ss*dbinom( 0:max(y) , prob=0.1 , size=10 )
points( c(0:max(y)+0.5) , p , pch=16 )

# show probability as a line
curve( 0.1 + 0.2*x , from=0 , to=3 , ylim=c(0,1.2) , ylab="probability" )
curve( 0.1 + 1*x , from=0 , to=3 , ylim=c(0,1) , ylab="probability" , add=TRUE )
lines( c(-1,4) , c(1,1) , lty=2 )
rect( -1 , 1 , 4 , 2 , col="red" , density=10 , angle=-45 )

# show probability as a logit
curve( 1/(1+exp(1 - 0.2*x)) , from=0 , to=3 , ylim=c(0,1.2) , ylab="probability" )
curve( 1/(1+exp(1 - 1*x)) , from=0 , to=3 , ylim=c(0,1) , ylab="probability" , add=TRUE )
curve( 1/(1+exp(1 - 2*x)) , from=0 , to=3 , ylim=c(0,1) , ylab="probability" , add=TRUE )
curve( 1/(1+exp(1 - 6*x)) , from=0 , to=3 , ylim=c(0,1) , ylab="probability" , add=TRUE )
rect( -1 , 1 , 4 , 2 , col="red" , density=10 , angle=-45 )
lines( c(-1,4) , c(1,1) , lty=2 )

# poisson examples

y <- rpois( 1000 , lambda=10 )
barplot( prep.bar(y) , names.arg=0:max(y) , xlab="y" , ylab="Frequency" )

y <- rpois( 30 , lambda=10 )
barplot( prep.bar(y) , names.arg=0:max(y) , xlab="y" , ylab="Frequency" )


# eagles example
library(MASS)
data(eagles)


# elephant example
d <- read.csv("case2201.csv")

# salamander example
d <- read.csv("case2202.csv")
barplot( prep.bar(d$SALAMAN) , names.arg=0:max(d$SALAMAN) , xlab="salamanders" , ylab="Frequency" )

m <- mle2( d$SALAMAN ~ dpois( a ) , start=list(a=3) )