# 1NRCAR on single subject's perio data:
# Code for setting up the Q matrix

# ordering of sites:  
# max then mand; ling then buc; tooth by tooth sort;  site by site sort

#Number of neighbors for each region:

num <- rep(4,168)
num[2+3*(0:55)] <- 2
num[c(1,42,43,84,85,126,127,168)] <- 2

#Define the neighbors:

adj <- matrix(NA,168,4)

# direct sites

l <- 2 + 3*(0:55) # 56 direct sites
for(i in 1:56){adj[l[i],1:2] <- c(l[i]-1,l[i]+1)}

# interprox sites

# left of interprox
l <- 3 + 3*(0:55) # there are 56, though some are end sites, 
                  # which are fixed below
#lingual
for(i in 1:14){adj[l[i],] <- c(l[i]-1,l[i]+1,l[i]+42,l[i]+43)}
for(i in 29:42){adj[l[i],] <- c(l[i]-1,l[i]+1,l[i]+42,l[i]+43)}
#buccal
for(i in 15:28){adj[l[i],] <- c(l[i]-1,l[i]+1,l[i]-42,l[i]-41)}
for(i in 43:56){adj[l[i],] <- c(l[i]-1,l[i]+1,l[i]-42,l[i]-41)}

# right of interprox
l <- 1 + 3*(0:55) # there are 56, though some are end sites, 
                  # which are fixed below
#lingual
for(i in 1:14){adj[l[i],] <- c(l[i]-1,l[i]+1,l[i]+41,l[i]+42)}
for(i in 29:42){adj[l[i],] <- c(l[i]-1,l[i]+1,l[i]+41,l[i]+42)}
#buccal
for(i in 15:28){adj[l[i],] <- c(l[i]-1,l[i]+1,l[i]-43,l[i]-42)}
for(i in 43:56){adj[l[i],] <- c(l[i]-1,l[i]+1,l[i]-43,l[i]-42)}

# end sites -- fixing the assignments just made
adj[1 ,]  <- c(2 ,43,NA,NA)
adj[42,]  <- c(41,84,NA,NA)
adj[43,]  <- c(1 ,44,NA,NA)
adj[84,]  <- c(83,42,NA,NA)
adj[85,]  <- adj[1 ,] + 84
adj[126,] <- adj[42,] + 84
adj[127,] <- adj[43,] + 84
adj[168,] <- adj[84,] + 84

adj <- c(t(adj))
adj <- adj[!is.na(adj)]

# Setting up the Q matrix

N <- 168
Q<-diag(num); count<-1
for(i in 1:N){for(j in 1:num[i]){Q[i,adj[count]]<--1; count<-count+1}}