library(RSiena)

# set object names
objnames <- paste0('saom.gofp.rate', rates, '.', saom.models)
subsets <- paste0('[\'mixed triad census\',]')
newobj <- paste0('gof', 1:length(saom.models))

for (i in 1:length(newobj)) eval(parse(text=paste0(newobj[i], ' <- ', objnames[i], subsets)))
rm(i, objnames, subsets, newobj)

# t-test table

saom.fit.ttests <- matrix(NA, 10, 8)
colnames(saom.fit.ttests) <- c('model1', 'model2', 'mean1', 'mean2', 'diff', 't-stat', 'p', 'df')

saom.fit.ttests[,'model1'] <- c(1, 2, 3, 4, 5, 6, 7, 6, 7, 7)
saom.fit.ttests[,'model2'] <- c(0, 1, 2, 3, 4, 4, 4, 5, 5, 6)

a <- t.test(gof1)
saom.fit.ttests[1,'mean1'] <- mean(gof1)
saom.fit.ttests[1,'mean2'] <- 0
saom.fit.ttests[1,'diff'] <- a$estimate    # diff
saom.fit.ttests[1,'t-stat'] <- a$statistic # t
saom.fit.ttests[1,'p'] <- a$p.value        # p
saom.fit.ttests[1,'df'] <- a$parameter     # df

for (i in 2:nrow(saom.fit.ttests)) {
  a <- eval(parse(text=paste0('t.test(gof', saom.fit.ttests[i,'model1'], 
                                   ', gof', saom.fit.ttests[i,'model2'], ', \'greater\', paired=T)')))
  saom.fit.ttests[i,'mean1'] <- eval(parse(text=paste0('mean(gof', saom.fit.ttests[i,'model1'], ')')))
  saom.fit.ttests[i,'mean2'] <- eval(parse(text=paste0('mean(gof', saom.fit.ttests[i,'model2'], ')')))
  saom.fit.ttests[i,'diff'] <- a$estimate    # diff
  saom.fit.ttests[i,'t-stat'] <- a$statistic # t
  saom.fit.ttests[i,'p'] <- a$p.value        # p
  saom.fit.ttests[i,'df'] <- a$parameter     # df
}
rm(i,a)

# print output to files if requested
if (print.results) {
  write.csv(round(saom.fit.ttests, 3), file='saom_fit_ttests.csv')
  write.csv(round(rbind(gof1, gof2, gof3, gof4, gof5, gof6, gof7), 4), 
            file='saom_fit_comparison.csv')
  write.csv(round(c(mean(gof1), mean(gof2), mean(gof3), mean(gof4), 
                    mean(gof5), mean(gof6), mean(gof7)), 4), 
            file='saom_fit_comparison_means.csv')
}

rm(gof1, gof2, gof3, gof4, gof5, gof6, gof7)