Description:

     'tables' calculates an estimate for the number of contingency
     tables with the counts satisfying a set of linear constraints. To
     implement the approximation, copy the files "tables.R" and
     "glm_cumulative_mod4.so" from the same directory as this file, to
     your working R directory.

Usage:

     tables(y, dimension, model, version, glm="geometric")

Arguments:

    y: a data vector containing all the counts of a contingency table.
   
    dimension: a vector with dimensions of the table.

    model:
          two-way table   :  Four models available: "I", "UA", "QI", "D".
          three-way table :  "(i)", "(ii)", "(iii)", "(iv)".
          
          (i)   corresponds to (R,C,Z)
          (ii)  corresponds to (RC,Z)
          (iii) corresponds to (RC,CZ)
          (iv)  corresponds to (RC,RZ,CZ).
                   

   version: 
         two-way table :   "Row", "Column"
         three-way table :  "IxJK", "JxIK", "KxIJ"
		 
   glm: 
       "geometric" - the default option for contingency tables with nonzero entries.
       "bernoulli" - for 0-1 tables.

Output:
$N0  - the estimate based on the double-saddlepoint approximation.
$N1  - the estimate based on the double-saddlepoint approximation 
       with linear correction term. 
$N2  - the estimate based on the double-saddlepoint approximation 
       with expontential correction term (recommended).
$ptm - elapsed time.

Details: see the reference paper.

References:
V.V.Zipunnikov, J.G.Booth and Ruriko Yoshida (2007), 
"Counting Tables using the Double Saddlepoint Approximation"

Example #1: sexual fun data

source("tables.R")
y<-c(7,2,1,2,7,8,5,8,2,3,4,9,3,7,9,14)
fit1.1<-tables(y,c(4,4,1),"I","Row")
fit1.2<-tables(y,c(4,4,1),"I","Column")

source("tables.R")
y<-c(7,2,1,2,7,8,5,8,2,3,4,9,3,7,9,14)
fit2.1<-tables(y,c(4,4,1),"UA","Row")
fit2.2<-tables(y,c(4,4,1),"UA","Column")

y<-c(7,2,1,2,7,8,5,8,2,3,4,9,3,7,9,14)
fit3.1<-tables(y,c(4,4,1),"QI","Row")
fit3.2<-tables(y,c(4,4,1),"QI","Column")

y<-c(7,2,1,2,7,8,5,8,2,3,4,9,3,7,9,14)
fit4.1<-tables(y,c(4,4,1),"D","Row")
fit4.2<-tables(y,c(4,4,1),"D","Column")

Example #2: pathologists' data

source("tables.R")
y<-c(22,5,0,0,0,2,7,2,1,0,2,14,36,14,3,0,0,0,7,0,0,0,0,0,3)
fit1.1<-tables(y,c(5,5,1),"I","Row")
fit1.2<-tables(y,c(5,5,1),"I","Column")

y<-c(22,5,0,0,0,2,7,2,1,0,2,14,36,14,3,0,0,0,7,0,0,0,0,0,3)
fit2.1<-tables(y,c(5,5,1),"UA","Row")
fit2.2<-tables(y,c(5,5,1),"UA","Column")

#Comment: Some cells are eliminated by conditioning in this example
y<-c(22,5,0,0,0,2,7,2,1,0,2,14,36,14,3)
fit3.1<-tables(y,c(5,3,1),"QI","Row")
fit3.2<-tables(y,c(5,3,1),"QI","Column")

y<-c(22,5,0,0,0,2,7,2,1,0,2,14,36,14,3,0,0,0,7,0,0,0,0,0,3)
fit4.1<-tables(y,c(5,5,1),"D","Row")
fit4.2<-tables(y,c(5,5,1),"D","Column")


Example #3: source("tables.R")

## Darwin Finch' data 
y<-c(
0,1,1,0,1,0,0,0,0,0,0,0,
0,1,1,0,1,0,0,0,0,0,0,0,
1,1,1,1,1,0,1,0,1,1,1,1,
1,1,1,1,0,0,1,0,1,1,1,1,
1,1,1,1,1,0,1,0,1,1,1,0,
1,1,1,0,1,0,1,0,1,1,0,0,
1,1,1,0,1,0,1,0,1,1,1,0,
1,1,1,1,1,0,1,0,1,1,1,0,
1,1,1,0,1,0,1,0,1,1,0,0,
1,1,1,1,1,0,0,0,1,1,1,0,
0,0,1,0,0,1,0,0,0,0,0,0,
1,1,1,1,1,0,1,1,1,1,0,0,
1,0,0,1,0,1,0,0,0,0,0,0,
1,1,1,0,1,0,1,0,0,1,0,0,
1,1,1,1,1,0,1,0,1,1,0,0,
1,0,0,1,0,0,0,0,0,0,0,0,
1,0,0,1,0,0,0,0,0,0,0,0)

fit1.1<-tables(y,c(12,17,1),"I","Row","bernoulli")
fit1.2<-tables(y,c(12,17,1),"I","Column","bernoulli")

Example#4: China cancer data

source("tables.R")
y<-c(126,35,100,61,908,497,688,807,913,336,747,598,235,58,172,121,
402,121,308,215,182,72,156,98,60,11,99,43,104,21,89,36)
fit1.1<-tables(y,c(2,2,8),"(i)","IxJK")
fit1.2<-tables(y,c(2,2,8),"(i)","JxIK")
fit1.3<-tables(y,c(2,2,8),"(i)","KxIJ")

y<-c(126,35,100,61,908,497,688,807,913,336,747,598,235,58,172,121,
402,121,308,215,182,72,156,98,60,11,99,43,104,21,89,36)
fit2.1<-tables(y,c(2,2,8),"(ii)","IxJK")
fit2.2<-tables(y,c(2,2,8),"(ii)","JxIK")
fit2.3<-tables(y,c(2,2,8),"(ii)","KxIJ")

y<-c(126,35,100,61,908,497,688,807,913,336,747,598,235,58,172,121,
402,121,308,215,182,72,156,98,60,11,99,43,104,21,89,36)
fit3.1<-tables(y,c(2,2,8),"(iii)","IxJK")
fit3.2<-tables(y,c(2,2,8),"(iii)","JxIK")
fit3.3<-tables(y,c(2,2,8),"(iii)","KxIJ")

y<-c(126,35,100,61,908,497,688,807,913,336,747,598,235,58,172,121,
402,121,308,215,182,72,156,98,60,11,99,43,104,21,89,36)
fit4.1<-tables(y,c(2,2,8),"(iv)","IxJK")
fit4.2<-tables(y,c(2,2,8),"(iv)","JxIK")
fit4.3<-tables(y,c(2,2,8),"(iv)","KxIJ")