**Wiekvoet**, and kindly contributed to R-bloggers)

Two weeks ago I showed rain data from six stations in Netherlands years 1906 till now. Last week I showed that frequency of days with and without rain differed between December 1906-1915 and December 2003-2012. This week I am considering the same data as t-distributed with a left truncation at 0 mm rain. This can be modeled very easy in JAGS. This model gives the posterior distribution of the amount of rain moved 1 mm of rain upwards and half a mm narrower, but both these intervals include 0.

### The data

### The model

The model is almost off the shelf. Truncation, two t-distributions with different means and standard deviations. Priors there taken from *Data Analysis Using Regression and Multilevel/Hierarchical Models*, Gelman and Hill.

sel1 <- all[(all$year <1916 | all$year>2003)

& all$Month==’Jan’

& all$location==’DE-BILT’

& all$day %in% seq(1,31,by=5),]

sel1$decade <- factor(c(‘1906-1915′,’2004-2013′)[1+(sel1$year>1950)])

datain <- list(

isCensored = 1-as.numeric(sel1$RD==0),

N = nrow(sel1),

RD = ifelse(sel1$RD==0,NA,sel1$RD/10),

period = 1+(sel1$year>1950)

)

model1 <- function() {

for ( i in 1:N ) {

isCensored[i] ~ dinterval( RD[i] , 0 )

RD[i] ~ dt( mu[period[i]] , tau[period[i]],nu )

}

for (i in 1:2) {

tau[i] <- 1/pow(sigma[i],2)

sigma[i] ~ dlnorm(sigmamu,sigmatau)

mu[i] ~ dnorm(mumu,mutau)

}

mumu ~ dnorm(0,1e-6)

mutau <- 1/pow(musigma,2)

musigma ~ dunif(0,100)

sigmamu ~ dnorm(0,1e-6)

sigmatau <- 1/pow(sigmasigma,2)

sigmasigma ~ dunif(1,100)

mudif <- mu[2]-mu[1]

sigmadif <- sigma[2]-sigma[1]

nu <- 1/nuinv

nuinv ~ dunif(0,.5)

}

params <- c(‘mu’,’sigma’,’mudif’,’sigmadif’,’nu’)

inits <- function() {

list(mumu = rnorm(1,0),

musigma = runif(1,1,2),

sigmamu = rnorm(1,0),

sigmasigma = runif(1,1,2),

nuinv = runif(1,0,.5))

}

jagsfit <- jags(datain, model=model1, inits=inits,

parameters=params,progress.bar=”gui”,

n.chains=1,n.iter=10000,n.burnin=3000,n.thin=3)

jagsfit

### Interpretation

**leave a comment**for the author, please follow the link and comment on their blog:

**Wiekvoet**.

R-bloggers.com offers

**daily e-mail updates**about R news and tutorials on topics such as: visualization (ggplot2, Boxplots, maps, animation), programming (RStudio, Sweave, LaTeX, SQL, Eclipse, git, hadoop, Web Scraping) statistics (regression, PCA, time series, trading) and more...