Revision 767

trunk/examples/admb-re/glmmadmb/glmmadmb.tpl (revision 767)
18 18
  init_ivector cor_block_start(1,M) 		// Not used: remove
19 19
  init_ivector cor_block_stop(1,M) 		// Not used: remove
20 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
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 22
  init_int link_type_flag   			// 0 log 1 logit 2 probit 3 inverse 4 cloglog 5 identity
23 23
  init_int rlinkflag                            // robust link function?
24 24
  init_int no_rand_flag   			// 0 have random effects 1 no random effects
25 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)
26 28
  // TESTING: remove eventually?
27 29
  init_int zi_kluge				// apply zi=0.001?
28 30
  init_int nbinom1_flag				// 1=NBinom1, 0=NBinom2
......
95 97

  
96 98
  // Determines "phases", i.e. when the various parameters  becomes active in the optimization process
97 99
  int pctr = 2;		// "Current phase" in the stagewise procedure
98
  // FIXME: move trunc_poisson earlier in like_type hierarchy (or add a scale parameter flag vector)
100
  // FIXME: move trunc_poisson, logistic earlier in like_type hierarchy (or add a scale parameter flag vector)
99 101
  int alpha_phase = like_type_flag>1 && like_type_flag!=6 ? pctr++ : -1;        // Phase 2 if active
100 102
  int zi_phase = zi_flag ? pctr++ : -1;                      			// After alpha
101 103
  int rand_phase = no_rand_flag==0 ? pctr++ : -1;    				// SD of RE's
......
412 414
      break;
413 415
    case 2:   // neg binomial
414 416
      if (cph<2)  // would like to use alpha_phase rather than 2 but it's in local_calcs
415
        tmpl = -square(log(1.0+y(_i,1))-log(1.0+lambda));
417
      	   tmpl = -square(log(1.0+y(_i,1))-log(1.0+lambda));
418
      if (cph<4)  
419
	  tmpl = -(1.0+y(_i,1))*square(log(1.0+y(_i,1))-log(1.0+lambda));
416 420
      else
417
        tmpl = log_negbinomial_density(y(_i,1),lambda,tau);
421
	  tmpl = log_negbinomial_density(y(_i,1),lambda,tau);
418 422
      break;
419 423
    case 3: // Gamma 
420 424
        tmpl = log_gamma_density(y(_i,1),alpha,alpha/lambda);
......
444 448
        tmpl = log_negbinomial_density(y(_i,1),lambda,tau)-log(1.0-pow(1.0+lambda/alpha,-alpha));
445 449
      break;
446 450
    case 8: // logistic 
447
      tmpl = -log(alpha) + (y(_i,1)-lambda)/alpha - 2*log(1+exp((y(_i,1)-lambda)/alpha));
451
	tmpl = -log(alpha) + (y(_i,1)-lambda)/alpha - 2*log(1+exp((y(_i,1)-lambda)/alpha));
448 452
      break;
449
    default:
453
    case 9: // beta-binomial
454
	Ni = sum(y(_i));
455
	tmpl = log_comb(Ni,y(_i,1)) + // log(C(Ni,y(_i,1)))
456
	    gammln(y(_i,2)+alpha*(1-lambda))+gammln(y(_i,1)+alpha*lambda)-gammln(Ni+alpha) + // lbeta(...)
457
	    -(gammln(alpha*(1-lambda))+gammln(alpha*lambda)-gammln(alpha)); // lbeta(...)
458
	break;
459
  default:
450 460
      cerr << "Illegal value for like_type_flag" << endl;
451 461
      ad_exit(1);
452 462
  }
......
454 464
  // Zero inflation part 
455 465
  // zi_kluge: apply ZI whether or not zi_flag is true
456 466
  if(zi_flag || zi_kluge) {
467
    // if (zi_model_flag) {
468
    // 
469
    // }
457 470
    if(y(_i,1)==0)
458 471
      g -= log(e2+pz+(1.0-pz)*mfexp(tmpl));
459 472
    else

Also available in: Unified diff