Revision 1027

trunk/examples/admb-re/glmmadmb/glmmadmb.tpl (revision 1027)
1
// 2011-03-24 version from Dave Fournier
2
// modified 2011-04-27 Hans Skaug, BMB
3
DATA_SECTION
4

  
5
  init_int n					// Number of observations
6
  init_int p_y					// Dimension of y(i) (multivariate response)
7
  init_matrix y(1,n,1,p_y)			// Observation matrix
8
  init_int p					// Number of fixed effects
9
  init_matrix X(1,n,1,p)			// Design matrix for fixed effects
10
  init_int M					// Number of RE blocks (crossed terms)
11
  init_ivector q(1,M)				// Number of levels of the grouping variable per RE block; Can be skipped
12
  init_ivector m(1,M)				// Number of random effects parameters within each block
13
  int sum_mq					// sum(m*q), calculated below: should be read from R
14
  init_int ncolZ				
15
  init_matrix Z(1,n,1,ncolZ)			// Design matrix for random effects
16
  init_imatrix I(1,n,1,ncolZ)			// Index vectors into joint RE vector "u" for each
17
  init_ivector cor_flag(1,M)			// Indicator for whether each RE block should be correlated
18
  init_ivector cor_block_start(1,M) 		// Not used: remove
19
  init_ivector cor_block_stop(1,M) 		// Not used: remove
20
  init_int numb_cor_params			// Total number of correlation parameters to be estimated
21
  init_int like_type_flag   			// 0 poisson 1 binomial 2 negative binomial 3 Gamma 4 beta 5 gaussian 6 truncated poisson 7 trunc NB 8 logistic 9 betabinomial
22
  init_int link_type_flag   			// 0 log 1 logit 2 probit 3 inverse 4 cloglog 5 identity
23
  init_int rlinkflag                            // robust link function?
24
  init_int no_rand_flag   			// 0 have random effects 1 no random effects
25
  init_int zi_flag				// Zero inflation (zi) flag: 1=zi, 0=no zi
26
  // init_int zi_model_flag                        // ZI varies among groups/covariates?
27
  // init_matrix G(1,n,1,ncolG)                    // Design matrix for zero-inflation (fixed effects)
28
  // TESTING: remove eventually?
29
  init_int zi_kluge				// apply zi=0.001?
30
  init_int poisshack                            // add e3 to poiss prob?
31
  init_int nbinom1_flag				// 1=NBinom1, 0=NBinom2
32
  init_int intermediate_maxfn			// Not used
33
  init_int has_offset				// Offset in linear predictor: 0=no offset, 1=with offset
34
  init_vector offset(1,n)			// Offset vector
35

  
36

  
37
  // Makes design matrix X orthogonal to improve numeric stability
38
  matrix rr(1,n,1,6)
39
  matrix phi(1,p,1,p)
40
  number ymax           // maximum y for bounding mean in nb and pois
41
 LOC_CALCS
42
  int i,j;
43
  phi.initialize();
44
  ymax=log(15.0*max(y)+1);
45
  for (i=1;i<=p;i++)
46
  {
47
    phi(i,i)=1.0;
48
  }
49
  dmatrix trr=trans(rr);
50
  trr(6).fill_seqadd(1,1);
51
  rr=trans(trr);
52

  
53
  dmatrix TX(1,p,1,n);
54
  TX=trans(X);
55
  for (i=1;i<=p;i++)
56
  {
57
    double tmp=norm(TX(i));
58
    TX(i)/=tmp;
59
    phi(i)/=tmp;
60
    for (j=i+1;j<=p;j++)
61
    {
62
      double a=TX(j)*TX(i);
63
      TX(j)-=a*TX(i);
64
      phi(j)-=a*phi(i);
65
    }
66
  }
67
  X=trans(TX);
68

  
69
  sum_mq = 0;
70
  for (i=1;i<=M;i++)
71
    sum_mq += m(i)*q(i);
72

  
73
  ofstream ofs("phi.rep");
74
  for (i=1; i<=p; i++)
75
  {
76
    for (j=1; j<=p; j++)
77
    {
78
      ofs << phi(i,j) << " ";
79
    }
80
    ofs << endl;
81
  }
82
  ofs << endl;
83

  
84
INITIALIZATION_SECTION
85
 tmpL 1.0
86
 tmpL1 0.0
87
 log_alpha 1
88
 pz .001
89

  
90
PARAMETER_SECTION
91
 LOC_CALCS
92

  
93
  // BMB: FIXME: do we need this?  formerly disallowed for binomial (was like_type_flag 2);
94
  //     would be problematic for binary data but otherwise OK.  Should test in R code
95
  //  if(zi_flag && (like_type_flag>=2))
96
  // {
... This diff was truncated because it exceeds the maximum size that can be displayed.

Also available in: Unified diff