*! VERSION 3.2.0, rewritten 6/14/92 for STATA version 3.0
/* Written by:
   D.H. Judson
   1013 James St.
   Newberg, OR 97132
   (503) 537-0870
   Program code enhancements and suggestions are welcome at the above address.
   This program calls ado files:
         stats.ado (STB-4, page 18)
         poisson.ado
   STATA version: 3.0 */
program define loglin
version 3.0
local options "Level(string) IRr"

if "`1'"=="" | substr("`1'",1,1)=="," {
	if "$S_E_cmd"!="poisson" { 
		error 301
	}
	parse "`*'"

	if "`level'"=="" {
		local level "95"
	}
	poisson, level(`level') `irr'
}

else {
local options /*
*/ "`options' Ltol(string) Iter(string) Offset(string) Fit(string) Anova Keep Resid COLLAPSE *"
local varlist "req ex"
local in "opt"
local if "opt"
local weight "opt fweight noprefix"
parse "`*'"
parse "`varlist'", parse(" ")

/* ASSOCIATE 1ST VAR==A, 2ND VAR==B, 3RD VAR==C, 4TH VAR==D */
mac def dv "`1'"
mac shift
mac def A "`1'"
mac def B "`2'"
mac def C "`3'"
mac def D "`4'"
mac def it=1
while ("${_$it}"~="") {
       mac def ${_$it}=substr("ABCD",$it,1)
       mac def it=$it+1
}

/* ERROR CHECK: MAXIMUM OF FOUR WAY TABLES ONLY */
cap assert "`5'"==""
if (_rc>0) {
      di in bl "ERROR: Only four independent variables can be specified.  You specified more than four."
      exit 9
}

/* DEFINE TOLERANCE */
if ("`ltol'"=="") {
      mac def _ltol "1e-7"
}

/* DEFINE ITERATIONS */
if ("`iter'"=="") {
      mac def _iter "100"
}

/* DEFINE OFFSET */
if ("`offset'"~="") {
      local poffset "+$_offset"
      local moffset "-$_offset"
      mac def offopt "offset(`offset')"
}

/* DEFINE CONFIDENCE INTERVAL */
if "`level'"=="" {
      local level "95"
}

/* DEFINE CELL WEIGHTS */
if ("`exp'"~="") {
      local wt "[fweight=`exp']"
}

/* SET DEFAULT MARGIN TO FIT IF NO MARGINS SPECIFIED */
if ("`fit'"=="") {
      di in bl "Note: Only the grand mean will be fit (one poisson parameter"
      di in bl "for all cells)."
      /* if no margins specified, fit grand mean only
         (i.e. one constant poisson parameter lambda
          for all levels/combinations of variables) */
}

/* CHECK COUNT VARIABLE */
cap assert ($dv>=0) & ($dv==int($dv)) `if' `in'
if (_rc>0) {
         di in bl "Note: You are responsible for interpretation of non-count dependent variable."
}

/* CODE TO CATCH -BREAK-, ERRORS, ETC., SEE P.305, RELEASE 3 REF. MANUAL */
if _N==0 {
      error 2000
}

local sfn "$S_FN"
tempfile user
quietly save "`user'"
capture {

/* TEMPORARILY COLLAPSE DATA SET USING STATS.ADO IF OPTION SPECIFIED */
if ("$_collaps"~="") {
      stats $dv `if' `in', sum by($A $B $C $D) nodescr collapse nowarning keep($dv)
      qui replace $dv=SM$dv `if' `in'
      qui drop SM$dv
}

/*  GENERATE ALL 1ST ORDER INDICATOR VARIABLES/VECTORS */
mac def i=1
while ("$_1"~="") {
      mac def newvar=substr("ABCD",$i,1)
      qui tab $_1, gen($newvar)

      if ("$_anova"~="") {          
         mac def k=_result(2) /* COUNT OF NUMBER OF CATEGORIES */
           mac def j=1
         mac def _retc=0
         while ($_retc==0) {
            qui cap replace $newvar$j=$newvar$j-$newvar$k
            mac def j=$j+1
            mac def _retc=_rc
                           /* ANOVA-LIKE CONSTRAINTS */
         }
         qui drop $newvar$k
      }
      else {
         qui drop ${newvar}1     /* REG-LIKE CONSTRAINTS--DEFAULT */
      }
      mac def i=$i+1
      macro shift
}

/* PARSE MARGINS TO BE FIT, STORE IN __$IT */ 
parse "`fit'", parse(",")
mac def it=1
while ("${_$it}"~="") {
      mac def __$it "${_$it}" 
      mac def it=$it+1
}

/* PARSE INDIVIDUAL INTERACTIONS, GENERATE ALL PARAMETERS TO FIT MARGINS */
mac def it=1
while ("${__$it}"~="") {
  if ("${__$it}"~=",") {
parse "${__$it}", parse(" ")
mac define z=("`1'"~="")+("`2'"~="")+("`3'"~="")+("`4'"~="")
                        /* number of configs to be fit */
mac define step=2               /* fit all 2-way configs first */
      while ($step<=4)*($step<=$z)*($z>1) {
         mac def i=1
         mac def j=2
         mac def k=3+($z<3)*2+($step<3)*2
         mac def l=4+($z<4)*2+($step<4)*2
         while ($i<=($z-$step+1))+($j<=($z-$step+2))+($k<=($z-$step+3))+($l<=($z-$step+4)) {
      
            /* DETERMINE LEVELS OF EACH MARGIN (VARIABLE) */
          qui cap tab ${$${_$i}} `if' `in'
          if _result(2)~=. {
            mac def M=(_rc==0)*(_result(2)-("`anova'"~="")) }
          else {
            mac def M=0 }
            
          qui cap tab ${$${_$j}} `if' `in'
          if _result(2)~=. {
            mac def N=(_rc==0)*(_result(2)-("`anova'"~="")) }
          else {
            mac def N=0 }

          qui cap tab ${$${_$k}} `if' `in'
          if _result(2)~=. {
            mac def O=(_rc==0)*(_result(2)-("`anova'"~="")) }
          else {
            mac def O=0 }

          qui cap tab ${$${_$l}} `if' `in'
          if _result(2)~=. {
            mac def P=(_rc==0)*(_result(2)-("`anova'"~="")) }
          else {
            mac def P=0 }

            /* GENERATE INDICATOR VARIABLES FOR INTERACTIONS */
            mac def m=0+($M>0)*(2-("`anova'"~="")) 
            mac def n=0+($N>0)*(2-("`anova'"~=""))
            mac def o=0+($O>0)*(2-("`anova'"~=""))
            mac def p=0+($P>0)*(2-("`anova'"~=""))
            mac def v1 "$${_$i}"
            mac def v2 "$${_$j}"
            mac def v3 "$${_$k}"
            mac def v4 "$${_$l}"

            while ($m<=$M)*($m~=0) {
              while ($n<=$N)*($n~=0) {
                  if ($o==0)*($p==0) {
                  cap gen byte $v1$v2$m$n=$v1$m*$v2$n
                  }
                  while ($o<=$O)*($o~=0) {
                  if ($p==0) {
                  cap gen byte $v1$v2$v3$m$n$o=$v1$m*$v2$n*$v3$o
                  }
                  while ($p<=$P)*($p~=0) {
                  cap gen byte $v1$v2$v3$v4$m$n$o$p=$v1$m*$v2$n*$v3$o*$v4$p
                  mac def p=0+($P>0)*($p+1)
                  }
                mac def o=0+($O>0)*($o+1)
                mac def p=0+($P>0)*(2-("`anova'"~=""))
                }
              mac def n=0+($N>0)*($n+1)
              mac def o=0+($O>0)*(2-("`anova'"~=""))
              }
            mac def m=0+($M>0)*($m+1)
            mac def n=0+($N>0)*(2-("`anova'"~=""))
            }
            /* END OF INTERACTION VARIABLES */ 
      
                 if $l<($z-$step+4) {
                          mac def l=$l+1 }
                  else if $k<($z-$step+3) {
                          mac def k=$k+1
                          mac def l=($l)*($step<4)+($k+1)*($step>=4) }
                 else if $j<($z-$step+2) {
                          mac def j=$j+1
                          mac def k=($k)*($step<3)+($j+1)*($step>=3) }
                 else if $i<($z-$step+1) {
                          mac def i=$i+1
                          mac def j=($j)*($step<2)+($i+1)*($step>=2) }
                 else mac def step=$step+1
         }
     }
  }
mac def it=$it+1
}
mac def vars " "
if "$A"~="" {noi di in gr "Variable $A =" in yel " A"
mac def vars "$vars A*"
}
if "$B"~="" {noi di in gr "Variable $B =" in yel " B"
mac def vars "$vars B*"
}
if "$C"~="" {noi di in gr "Variable $C =" in yel " C"
mac def vars "$vars C*"
}
if "$D"~="" {noi di in gr "Variable $D =" in yel " D"
mac def vars "$vars D*"
}
if "`fit'"=="" {noi di in gr "Margins fit: Grand mean only"
}
else noi di in gr "Margins fit: " in yel "`fit'"

if "`anova'"~="" {
   noi di in gr "Note: Anova-like constraints are assumed.  The last level of each"
   noi di in gr "variable (and all interactions with it) will be dropped from estimation."
   noi di in gr "The variable codings are constrained to sum to zero, so the last"
   noi di in gr "level will equal -1 times the sum of the other levels."
}
else {
   noi di in gr "Note: Regression-like constraints are assumed.  The first level of each"
   noi di in gr "variable (and all iteractions with it) will be dropped from estimation."
}

if "`fit'"=="" {
         noi poisson $dv `in' `if' `wt', ltol(`ltol') iter(`iter') $offopt level(`level') `irr'
}
else noi poisson $dv $vars `in' `if' `wt', ltol(`ltol') iter(`iter') $offopt level(`level') `irr'

if "`resid'"~="" {
      qui predict __lcnts `in' `if'
      gen cellhat=exp(__lcnts) `in' `if'
            cap replace cellhat=0 if (`exp'==0)|(`exp'==.)|($dv==.) `in'
      gen resid=$dv-cellhat `in' `if'
      gen stdres=resid/sqrt(cellhat) `in' `if'
      lab var resid "RESIDUAL"
      lab var stdres "STANDARDIZED RESIDUAL"
      lab var cellhat "ESTIMATED EXPECTED CELL FREQ."
      qui sort $A $B $C $D $dv
      format $_varlist %3.0f
      format cellhat %8.3f
      format resid %7.3f
      format stdres %7.3f
      noi list `varlist' cellhat resid stdres, nodisplay noobs 
      format `varlist' %9.0g
      qui drop __lcnts
}

if "`keep'"=="" {
      qui drop $vars
}
else if "`collaps'"~=""{
      noi di in bl "NOTE:  Cannot keep your results.  The table is collapsed."
}
else {
      save "`user'", replace
}
}
local rc=_rc
qui cap use "`user'", clear
mac def S_FN "`sfn'"
cap erase "`user'"
error `rc'
}
end