Hi All,
I'm fitting a kind of VAR model (vector autoregression), based on a
Lotka-Volterra scheme, to population time series data for a set of "i"
species. The model works pretty fine, but now I want to estimate the
eigenvalues of the Jacobian of the interaction matrix, which I guess is
analogous to the matrix of AR coefficients in a standard VAR model.
Does anyone has an example on how to do this with the eigen.vals(X)
function in OpenBUGS? I guess I have to assemble a square matrix "X",
but I'm really doing something wrong.
Below I paste the code for the model, including the estimation of the
individual terms within the Jacobian matrix.
Thank you very much in advance for your kind help,
Regards,
Pablo Almaraz
#####################
### BEGIN CODE: #####
model {
# Lotka-Volterra model
for(t in 2:NYears){
for(i in 1:NSpecies){
n.mu[t,i] <- max(n[t-1,i] + r[i]*(1 -
inprod(alpha[i,],n[t-1,])/k[i]), 0)
}
n[t,1:NSpecies] ~ dmnorm(n.mu[t,1:NSpecies], Tau[,])
}
## Interaction coefficients
for(i in 1:NSpecies){
alpha[i,i] <- 1
for(j in i+1:NSpecies){
alpha[i,j] ~ dnorm(0, prec.tau.alpha[i,j])
prec.tau.alpha[i,j] <- 1/tau.alpha2[i,j]
tau.alpha2[i,j] <- pow(tau.alpha[i,j], 2)
tau.alpha[i,j] ~ dunif(0, 3)
alpha[j,i] ~ dnorm(0, prec.tau.alpha[j,i])
prec.tau.alpha[j,i] <- 1/tau.alpha2[j,i]
tau.alpha2[j,i] <- pow(tau.alpha[j,i], 2)
tau.alpha[j,i] ~ dunif(0, 3)
}
}
## Setup the entries for the Jacobian matrix:
for(i in 1:NSpecies){
J[i,i] <- 1-(r[i]/k[i])
for(j in i+1:NSpecies){
J[i,j] <- -(alpha[i,j]*r[i])/k[i]
J[j,i] <- -(alpha[j,i]*r[j])/k[j]
}
}
## Residual correlations
for(i in 1:NSpecies){
for(j in 1:NSpecies){
Corr.sigma[i,j] <- Sigma[i,j]/sqrt(Sigma[i,i]*Sigma[j,j])
}
n.sigma2[i] <- Sigma[i,i]
r[i] ~ dnorm(0, 0.1) I(0,)
k[i] ~ dnorm(0, 0.01) I(0,)
}
# Residual environmental covariance matrix (Sigma)
Tau[1 : NSpecies,1 : NSpecies] ~ dwish(Rho[ , ], NSpecies)
Sigma[1:NSpecies,1:NSpecies] <- inverse(Tau[,])
}
#####################
### END CODE: #######
-------------------------------------------------------------------
This list is for discussion of modelling issues and the BUGS software.
For help with crashes and error messages, first mail [log in to unmask]
To mail the BUGS list, mail to [log in to unmask]
Before mailing, please check the archive at www.jiscmail.ac.uk/lists/bugs.html
Please do not mail attachments to the list.
To leave the BUGS list, send LEAVE BUGS to [log in to unmask]
If this fails, mail [log in to unmask], NOT the whole list
|