rjmcmc: Model discriminiation in a Bayesian context for spatial NBDA.

In spatialnbda: Performs spatial NBDA in a Bayesian context

Description

An RJMCMC algorithm is used to achieve model discrimination between the null model which contains only the baseline rate parameter and the full model which contains both the baseline rate and social parameters.

Usage

rjmcmc(formatteddata, its, pilot_tuner1, pilot_tuner2, start1, start2, start3, p1, p2)

Arguments

formatteddata	Formatted data using the FormatData function.
its	Number of iterations
pilot_tuner1	Tuner for proposal distribution for the social parameter.
pilot_tuner2	Tuner for the proposal distribution for the baseline rate parameter.
start1	Start value for the social parameter
start2	Start value for the baseline rate parameter
start3	Start model
p1	Unifor prior variance tuner for the baseline rate
p2	Uniform prior variance tuner for the social parameter

Details

It is important to check that the chains have mixed which using this function. A rough way would be to view the trace plots printed.

Value

The output is a table with the number of iterations for which the Markov chain spent in each visited model.

Examples

#Example 1
data(timearray)
data(idarray)
data(socialx)
data(socialy)

Times = timearray[,1]
Ids = idarray[,1]
lenh = length(Times)
Groups = rep(1,length(Times))
Events = c(1:length(Times))

socialites = matrix(1,nrow=lenh,ncol=lenh)
x = socialx
y = socialy

plot(x[,1],y[,1],xlab="x",ylab="y",cex=2,pch=16,main="Point pattern of nest positions")

areas = calculate.areas(x[,1],y[,1],rep(0.2,lenh),1000)
spatialareas = areas
len = length(x[,1])
Diffusions = rep(1,len)
for(i in 2:10){
  addon = rep(i,len)
  Diffusions = c(Diffusions,addon)
  
}


spatialnetwork = matrix(0,nrow=len,ncol=len)
for(i in 1:len){
  for(j in i:len){ 
    template = spatialareas[[i]][j]
    spatialnetwork[i,j] = spatialnetwork[j,i] = template
    #spatialareas[[i]]=NULL
    
  }
  
}


shape = FormatData(Times,spatialnetwork,Ids,Groups,Diffusions,Events,spatialnetwork)



#ptm <- proc.time()
#mcmc(shape,10000,0.05,0.05,-3,-5)
#proc.time() - ptm

#ptm <- proc.time()
#rjmcmc(shape,10000,5,1,-3,-3,1,10,10)
#proc.time() - ptm


# Example 2
data(papertimes)
data(papernests)
data(x)
data(y)
z = array(0,c(length(x[,1]),1))# setting up array for storing spatial covariate information

for(i in 1:70){   # simulating spatial covariate information
xx = x[,1][i]
yy = y[,1][i]
z[i] = (3*xx + 14*yy) * exp(2 * (.4*xx - 1)) 
}



Times = papertimes[,1]
Ids = papernests[,1]
Diffusions = rep(1,length(Times))
Groups = rep(1,length(Times))
Events = c(1:length(Times))
socialites = matrix(1,nrow=70,ncol=70)


plot(x[,1],y[,1],xlab="x",ylab="y",cex=2,pch=16,main="Point pattern of nest positions")



areas = calculate.areas(x[,1],y[,1],rep(0.05,70),1000)
spatialareas = areas
len = length(x[,1])

spatialnetwork = matrix(0,nrow=len,ncol=len)
for(i in 1:len){
  for(j in i:len){ 
    template = spatialareas[[i]][j]
    spatialnetwork[i,j] = spatialnetwork[j,i] = template
    #spatialareas[[i]]=NULL
    
  }
  
}


shape = FormatData(Times,spatialnetwork,Ids,Groups,Diffusions,Events,spatialnetwork,z)


#ptm <- proc.time()
#mcmc(shape,10000,5,1,-5,-6)
#proc.time() - ptm

#ptm <- proc.time()
#nullmcmc(shape,10000,1,-5)
#proc.time() - ptm


#ptm <- proc.time()
#rjmcmc(shape,10000,5,1,0,0,2,5,5)
#proc.time() - ptm

spatialnbda documentation built on May 2, 2019, 8:54 a.m.