Skip to content

Instantly share code, notes, and snippets.

@mortonjt
Last active November 26, 2019 00:30
Show Gist options
  • Save mortonjt/227a0058194c2fa0ecf06807e1315d35 to your computer and use it in GitHub Desktop.
Save mortonjt/227a0058194c2fa0ecf06807e1315d35 to your computer and use it in GitHub Desktop.
This performs a very simple negative binomial regression tailored for differential abundance analysis
data {
int<lower=0> N; // number of samples
int<lower=0> D; // number of dimensions
int<lower=0> p; // number of covariates
real depth[N]; // sequencing depths of microbes
matrix[N, p] x; // covariate matrix
int y[N, D]; // observed microbe abundances
}
parameters {
// parameters required for linear regression on the species means
matrix[p, D-1] beta;
real reciprocal_phi;
}
transformed parameters {
matrix[N, D-1] lam;
matrix[N, D] lam_clr;
matrix[N, D] prob;
vector[N] z;
real phi;
phi = 1. / reciprocal_phi;
z = to_vector(rep_array(0, N));
lam = x * beta;
lam_clr = append_col(z, lam);
}
model {
// setting priors ...
reciprocal_phi ~ cauchy(0., 5.);
for (i in 1:D-1){
for (j in 1:p){
beta[j, i] ~ normal(0., 5.); // uninformed prior
}
}
// generating counts
for (n in 1:N){
for (i in 1:D){
target += neg_binomial_2_log_lpmf(y[n, i] | depth[n] + lam_clr[n, i], phi);
}
}
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment