> d
   dept applicant_gender admit reject applications male dept_id
1     A             male   512    313          825    1       1
2     A           female    89     19          108    0       1
3     B             male   353    207          560    1       2
4     B           female    17      8           25    0       2
5     C             male   120    205          325    1       3
6     C           female   202    391          593    0       3
7     D             male   138    279          417    1       4
8     D           female   131    244          375    0       4
9     E             male    53    138          191    1       5
10    E           female    94    299          393    0       5
11    F             male    22    351          373    1       6
12    F           female    24    317          341    0       6

----------chapter 10-----------------

> precis(m10.6)   
    mean   sd  5.5% 94.5% n_eff Rhat
a  -0.82 0.05 -0.90 -0.74    92    1
bm  0.60 0.07  0.49  0.70    94    1

> stancode(m10.6)
data{
    int applications[12];
    int admit[12];
    real male[12];
}
parameters{
    real a;
    real bm;
}
model{
    vector[12] p;
    bm ~ normal( 0 , 10 );
    a ~ normal( 0 , 10 );
    for ( i in 1:12 ) {
        p[i] = a + bm * male[i];
        p[i] = inv_logit(p[i]);
    }
    admit ~ binomial( applications , p );
}

# Model in logpdf form (written by Scott Spencer)

m13.2_model = "
  data {
    int N;
    int<lower=0, upper=1> male[N];
    int N_admit[N];
    int N_applications[N];
  }
  parameters {
    real a;
    real bm;
  }
  model {
    vector[N] p_admit;
    for (i in 1:N) p_admit[i] = a + bm * male[i];
    target += binomial_logit_lpmf(N_admit | N_applications, p_admit);
    target += cauchy_lpdf(a | 0, 5);
    target += normal_lpdf(bm | 0, 5);
  }
";

----------chapter 13-----------------

> precis( m13.2 , depth=2 ) # depth=2 to display vector parameters
            mean   sd  5.5% 94.5% n_eff Rhat
a_dept[1]   0.68 0.10  0.51  0.83  7886    1
a_dept[2]   0.63 0.11  0.45  0.81  8279    1
a_dept[3]  -0.58 0.08 -0.70 -0.46 10468    1
a_dept[4]  -0.62 0.09 -0.76 -0.48  8972    1
a_dept[5]  -1.06 0.10 -1.22 -0.91 12411    1
a_dept[6]  -2.61 0.16 -2.87 -2.37 10834    1
a          -0.58 0.65 -1.57  0.41  6458    1
bm         -0.10 0.08 -0.22  0.03  6646    1
sigma_dept  1.48 0.58  0.85  2.48  5787    1

> stancode(m13.2)
data{
    int applications[12];
    int admit[12];
    real male[12];
    int dept_id[12];
}
parameters{
    vector[6] a_dept;
    real a;
    real bm;
    real<lower=0> sigma_dept;
}
model{
    vector[12] p;
    sigma_dept ~ cauchy( 0 , 2 );
    bm ~ normal( 0 , 1 );
    a ~ normal( 0 , 10 );
    a_dept ~ normal( a , sigma_dept );
    for ( i in 1:12 ) {
        p[i] = a_dept[dept_id[i]] + bm * male[i];
        p[i] = inv_logit(p[i]);
    }
    admit ~ binomial( applications , p );
}

# Varying intercept logpdf model (written by Scott Spencer)

m13.2_model = "
  data {
    int N;
    int n[N];
    int m[N];
    int A[N];
    int N_depts;
    int dept[N];
  }
  parameters {
    real alpha;
    real beta;
    real<lower=0> sigma;
    vector[N_depts] a_dept;
  }
  model {
    vector[N] p;
    target += normal_lpdf(a_dept | alpha, sigma);
    target += normal_lpdf(alpha | 0, 10);
    target += normal_lpdf(beta | 0, 1);
    target += cauchy_lpdf(sigma | 0, 2);
    for(i in 1:N) p[i] = a_dept[dept[i]] + beta * m[i];
    target += binomial_logit_lpmf(A | n, p);
  }
  generated quantities {
    vector[N] log_lik;
    {
      vector[N] p;
      for (i in 1:N) {
        p[i] = a_dept[dept[i]] + beta * m[i];
        log_lik[i] = binomial_logit_lpmf(A[i] | n[i], p[i]);
      }
    }
  }
";
