library(RSiena)
source('00_functions.R')

# rename result objects from selected model specification
for (i in village.codes) {
  
  objname <- paste0('saom.', i, '.rate', rates, '.', saom.model)
  assign(paste0('saom.result.', i), eval(parse(text=objname)))
  
}
rm(i, objname)

# p-value tables
saom.gofp <- matrix(NA, 9, length(village.codes))
rownames(saom.gofp) <- c('talk: outdegree', 'talk: indegree', 'talk: triad census',
                         'talk: geodesics', 'talk: clique census', 'use: outdegree',
                         'use: indegree', 'mixed triad census', 'mixed triad census - agreement')

# gofs
for (i in 1:length(village.codes)) {
  eval(parse(text=paste0('s.result <- saom.result.v', i)))
  saom.gofp[1,i] <- sienaGOF(s.result, OutdegreeDistribution, 
                             varName='talk', levls=0:8, tested=F)$Joint$p
  saom.gofp[2,i] <- sienaGOF(s.result, IndegreeDistribution, 
                             varName='talk', levls=0:11, tested=F)$Joint$p
  saom.gofp[3,i] <- sienaGOF(s.result, TriadCensus, 
                             varName='talk', tested=F)$Joint$p
  saom.gofp[4,i] <- sienaGOF(s.result, GeodesicDistribution, 
                             varName='talk', tested=F)$Joint$p
  saom.gofp[5,i] <- sienaGOF(s.result, CliqueCensus, 
                             varName='talk', tested=F)$Joint$p
  saom.gofp[6,i] <- sienaGOF(s.result, OutdegreeDistribution, 
                             varName='use', levls=0:7, tested=F)$Joint$p
  saom.gofp[7,i] <- sienaGOF(s.result, IndegreeDistribution, 
                             varName='use', levls=0:40, tested=F)$Joint$p
  saom.gofp[8,i] <- sienaGOF(s.result, mixedTriadCensus, 
                             varName=c('talk', 'use'), tested=F)$Joint$p
  saom.gofp[9,i] <- sienaGOF(s.result, mixedTriadCensus.agree, 
                             varName=c('talk', 'use'), tested=F)$Joint$p
  print(paste0(saom.model, ', village ', i, ': done!'))
}
rm(i, s.result)

# rename gof p-values object to reflect rate parameters and specification
objname <- paste0('saom.gofp.rate', rates, '.', saom.model)
assign(objname, saom.gofp)

# print results to file if requested
if (print.results) write.table(saom.gofp, file=paste0(objname, '.csv'), sep=',')

# remove unnecesary objects
rm(saom.gofp, objname, saom.model)
rm(list=ls(pattern='saom.result'))
rm(siena07ToConvergence, household.aggregate, predictive.acc, reg.backward,
   GeodesicDistribution, CliqueCensus, mixedTriadCensus.agree)