I tried to convert your model to my case (simplified).
I assume from your documentation and discourse reply:
- your targets == my genes
- your regulators == my covariates (I just have one in this example)
- if you have 1 regulator/covariate: relative_weights == [1, 1, …, 1], so we can ignore
- if the sign of the slope is not supposed to not cause non-identifiability issues, regulation_signs is omitted and sd_regulatory_input must exist in the domain -(Inf, Inf)
My basic link function with equation is
inv_logit( intercept + X * beta ) * k
Given all this my model is
data {
int<lower = 0> T; // tube colleciton of genes for a patient
vector<lower = 0>[T] y; // RNA-seq counts
vector[T] X; // design matrix
}
parameters {
real mean_regulatory_input;
real sd_regulatory_input;
real<lower=0> sigma;
real mean_observed_value;
}
transformed parameters {
real observed_sd_input = sd(X);
real observed_mean_input = mean(X);
real slope = sd_regulatory_input / observed_sd_input;
real intercept = mean_regulatory_input - slope * observed_mean_input;
real k = mean_observed_value / mean(inv_logit(intercept + X * slope));
}
model {
for (t in 1:T) y[t] ~ normal( inv_logit( intercept + X[t] * slope ) * k, sigma );
mean_regulatory_input ~ normal(0,5);
sd_regulatory_input ~ normal(0,5);
sigma ~ cauchy(0,2.5);
mean_observed_value ~ normal(0,5);
}
it does not converge so well although tends to predict reasonable values for parameters.
My questions are: it is correct in your opinion? Am I missing something on the prior usage?
Slope == 0
Slope == 5
