Skip to content

Instantly share code, notes, and snippets.

Embed
What would you like to do?
Stan model for noisy observations of a sigmoidal function, where the noise is truncated to a finite interval
functions {
vector ft(vector t, real r, real y0, real tC, real L) {
vector[num_elements(t)] mu;
vector[num_elements(t)] exprt;
exprt = exp(r*(t-tC));
mu = y0*L*exprt ./ (L + y0*(exprt - 1));
return mu;
}
vector dfdt(vector t, real r, real tC, real L) {
vector[num_elements(t)] dmu;
vector[num_elements(t)] sddenom;
for (i in 1:num_elements(t)) {
sddenom[i] = ((exp(r*tC) + exp(r*t[i]))^(-2));
}
dmu = r*L*exp(r*(t+tC)) .* sddenom;
return dmu;
}
}
data {
int<lower = 1> M;
int<lower = 1> N;
real<lower = 1> maxY;
real tcrit;
matrix<lower=0, upper=maxY>[M,N] y;
vector[M] t;
real eps;
int<lower = 1> P;
vector[P] tpred;
}
parameters {
real<lower = 0> r;
real<lower = 0, upper=maxY> mu0;
}
transformed parameters {
vector[M] curr_mu;
vector[M] curr_sd;
curr_mu = ft(t, r, mu0, tcrit, maxY);
curr_sd = dfdt(t, r, tcrit, maxY) + eps;
}
model {
r ~ gamma(1, 10);
mu0 ~ uniform(0, maxY);
for (i in 1:M) {
y[i,] ~ normal(curr_mu[i], curr_sd[i]);
}
}
generated quantities{
vector[P] pred_mu;
vector[P] pred_sd;
pred_mu = ft(tpred, r, mu0, tcrit, maxY);
pred_sd = dfdt(tpred, r, tcrit, maxY) + eps;
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
You can’t perform that action at this time.