model {
 	# Likelihood
	for (i in 1:Nareas) {
		for (k in 1:Ndiseases) {
			Y[i, k] ~ dpois(mu[i, k])
			log(mu[i, k]) <- log(E[i, k]) + alpha[k] + S[i,k]
  			RR[i,k] <- exp(alpha[k] + S[i,k])
		}
  	}

  	for (i in 1:Nareas){
		S[i,1]<-Het[1,i]*CholDis1[1,1]+Spatial[1,i]*CholDis2[1,1]
		for(k in 2:Ndiseases){
			S[i,k]<-inprod(Het[1:k,i],CholDis1[k,1:k])+inprod(Spatial[1:k,i],CholDis2[k,1:k])
		}
	}

	for(j in 1:Ndiseases){
		Spatial[j,1:Nareas]~car.normal(adj[], weights[], num[],1)
		for(i in 1:Nareas){Het[j,i]~dnorm(0,1)}
	}

	#Between-diseases structure
	#Cholesky triangle of a general correlation matrix
  	#diagonal
	CholDis1[1,1]<-sdHet[1]
	CholDis2[1,1]<-sdSpatial[1]
	for(j in 2:Ndiseases){
		diag1[j]<-pow(sdHet[j],2)-inprod(CholDis1[j,1:(j-1)],CholDis1[j,1:(j-1)])
		CholDis1[j,j]<-sqrt(abs(diag1[j]))
		diag2[j]<-pow(sdSpatial[j],2)-inprod(CholDis2[j,1:(j-1)],CholDis2[j,1:(j-1)])
		CholDis2[j,j]<-sqrt(abs(diag2[j]))
	}

	#Lower triangle
	#First column
	for(i in 2:Ndiseases){
		CholDis1[i,1]<-sdHet[i]*roDis[i,1]
		CholDis2[i,1]<-sdSpatial[i]*roDis[i,1]
	}
	#Rest of columns
	for(j in 2:(Ndiseases-1)){
		for(i in (j+1):Ndiseases){
			CholDis1[i,j]<-(sdHet[i]*sdHet[j]*roDis[i,j]-inprod(CholDis1[i,1:(j-1)],CholDis1[j,1:(j-1)]))/CholDis1[j,j]
			CholDis2[i,j]<-(sdSpatial[i]*sdSpatial[j]*roDis[i,j]-inprod(CholDis2[i,1:(j-1)],CholDis2[j,1:(j-1)]))/CholDis2[j,j]
		}
	}

	for(j in 1:Ndiseases){roDis[1,j]<-0}
	for(i in 2:Ndiseases){for(j in 1:(i-1)){roDis[i,j]~dunif(-1,1)}
	for(j in i:Ndiseases){roDis[i,j]<-0}}

	#Condition of being positive defined the covariance matrix
	uno<-1 
	uno~dbern(condition)
	condition<-step(sum(subcondition[2:Ndiseases])-(Ndiseases-1))
	for(j in 2:Ndiseases){subcondition[j]<-step(diag1[j])}

	# Other priors
	for (k in 1:Ndiseases) {
		alpha[k] ~ dflat()
		sdSpatial[k] ~ dunif(0,10)
		sdHet[k] ~ dunif(0,10)
	}
}