florianhartig

# Praktische Vorlesung, Einführung Statistik, WS 13/14
# Florian Hartig, http://www.biom.uni-freiburg.de/mitarbeiter/hartig

?read.csv

# Wiederholung likelihood

# Als estes erstellen wir einen leeren Plot

plot(NULL, NULL, xlim=c(-4,6), ylim = c(0,0.5), ylab = "Wahrscheinlichkeit(sdichte)", xlab = "Observed value")

florianhartig

library(IDPmisc)

panel.hist <- function(x, ...)
{
  usr <- par("usr"); on.exit(par(usr))
  par(usr = c(usr[1:2], 0, 1.5) )
  h <- hist(x, plot = FALSE)
  breaks <- h$breaks; nB <- length(breaks)
  y <- h$counts; y <- y/max(y)
  rect(breaks[-nB], 0, breaks[-1], y, col="blue4", ...)

florianhartig

# plot layout
x=12
y=12
# number of diff samples
n=3

# plot A, which holds
A=array(NA,c(x,y))
# repeat 1 through n as many times as subplots are in A
s=rep(1:n, x*y/n)

florianhartig

library(mvtnorm)   # to draw multivariate normal outcomes
library(R2jags)    # JAGS-R interface

# function that makes distance matrix for a side*side 2D array  
dist.matrix <- function(side)
{
  row.coords <- rep(1:side, times=side)
  col.coords <- rep(1:side, each=side)
  row.col <- data.frame(row.coords, col.coords)
  D <- dist(row.col, method="euclidean", diag=TRUE, upper=TRUE)

florianhartig

rm(list=ls(all=TRUE))

#Analysis of likelihood, p-value, and Bayes for binomial model,
#10 trials, 3 success, unknown coin, want to do inference

trials = 10
success = 3

# GETTING THE MLE ESTIMATE

florianhartig

rm(list=ls(all=TRUE))

library(rjags)

# assuming the data is created from an ecological system that creates an
# exponential size distribution (e.g. you sample individuals from a population that can be 
# expected to follow this distribution), but this measurments are done with
# a considerable lognormal observation error
# for a realistic application see http://www.plosone.org/article/info%3Adoi%2F10.1371%2Fjournal.pone.0058036

florianhartig

# modified from http://stackoverflow.com/questions/1401904/painless-way-to-install-a-new-version-of-r-on-windows

# run on old computer / r version
setwd("/Users/Florian/Dropbox/temp") # or any other existing temp directory
packages <- installed.packages()[,"Package"]
save(packages, file="Rpackages")


# run on new computer / r version
setwd("/Users/Florian/Dropbox/temp") # or any other existing temp directory

florianhartig

# from https://tonybreyal.wordpress.com/2011/11/24/source_https-sourcing-an-r-script-from-github/
# and http://stackoverflow.com/questions/7715723/sourcing-r-script-over-https

source_https <- function(url, ...) {
  # load package
  require(RCurl)
 
  # parse and evaluate each .R script
  sapply(c(url, ...), function(u) {
    eval(parse(text = getURL(u, followlocation = TRUE, cainfo = system.file("CurlSSL", "cacert.pem", package = "RCurl"))), envir = .GlobalEnv)

florianhartig

## Load libraries
library(mclogit)
library(reshape)
library(rjags)
library(R2WinBUGS) ## for write.model

## Load the data file from mclogit
data(Transport)

## Do the mclogit model

florianhartig

############ CREATE ZERO-INFLATED GLMM DATA #################

# This first part creates a dataset with beetles counts across an altitudinal gradient (several plots each observed several years), with a random intercept on year and zero-inflation. 

altitude = rep(seq(0,1,len = 50), each = 20)
dataID = 1:1000
spatialCoordinate = rep(seq(0,30, len = 50), each = 20)

# random effects + zeroinflation
plot = rep(1:50, each = 20)