[R-sig-ME] parameter restrictions

[dave fournier]

I don't think it is very difficult. Something like this should do the job.
It depends on how you want to model the correlation of the random effects.
Perhaps the ADMB list would be the place to discuss it.



DATA_SECTION
   int nfixed
  !! nfixed=7;  // for this example
   init_int nobs   // number of observations
   init_int nrand   // number of observations
   init_vector Y(1,nobs)   // observations
   init_matrix X(1,nobs,1,nfixed)  // fixed effects design matrix
   init_matrix Z(1,nobs,1,nrand)   // random effects design matrix
PARAMETER_SECTION
   init_number a
   init_bounded_vector alpha(2,nfixed,0.0,1.0)
   objective_function_value f
  !! int nr1=nrand*(nrand-1)/2; // number of parameters needed for 
correlation
   init_vector c(1,nr1)   // parameters for choleski decomp
                                        // of the correlation matrix
   init_vector log_sigma(1,nrand)           // log of variances
   random_effects_vector u(1,nrand)
PROCEDURE_SECTION
   // parameterize the b's to satisfy constraints
   dvar_vector b(1,nfixed);
   b(1)=-exp(a);
   dvariable ssum=0.0;
   for (int i=2;i<=nfixed;i++) {
     ssum-=b(i-1);
     b(i)=alpha(i)*ssum;
   }
   // choleski decomp of correlation matrix
   dvar_matrix C(1,nrand,1,nrand);
   C.initialize();
   int ioffset=1;
   for (int i=1;i<=nrand;i++) {
     C(i,i)=1;
     if (i>1)
     {
       C(i)(1,i-1)=c(ioffset,ioffset+i-1).shift(1);
       C(i)/=norm(C(i));
       ioffset+=i-1;
     }
     C(i)*=exp(log_sigma(i));   // turn it into Choleski decomp of 
covariance
                                // matrix
   }
   dvar_vector v=C*u;           // now v has right covariance
   dvar_vector pred_Y=X*b+Z*v;  // predicted Y
   dvar_vector r=Y-pred_Y;
   dvariable w=norm2(r)/nobs;   // mean square error
   f=0.5*nobs*log(w);    // concemntrated likelihood
   f+=0.5*norm2(u);    // log prior  for u