Hi,

This is my code:

evol_tree_code = {

‘data {’

’ int<lower=0> T; ‘

’ int<lower=0> N; ‘

’ real<lower=0,upper=1> Fhat[N,T]; ‘

’ matrix<lower=0,upper=1> [N,N] UTree; ‘

’}’

‘parameters {’

’ real<lower= 0,upper=0.3> mu[N]; ‘

’ vector<lower=0> [N] delta_tau; ‘

’}’

‘transformed parameters {’

’ matrix <lower=0> [N,T] numM; ‘

’ matrix <lower=0,upper=1> [N,T] M; ‘

’ matrix <lower=0> [N,T] F; ‘

’ vector<lower=0> [N] tau;’

’ real<lower=0> den[T];’

’ tau = UTree’‘*delta_tau; ‘*

‘for(t in 1:T){’

’ den[t] = 0.001;//this should be a small value’

’ for (n in 1:N){’

’ if (tau[n] > t - 1){’

’ numM[n,t] = 0.001;’

’ }else{’

’ numM[n,t] = exp(mu[n](t-1 - tau[n]));’

’ }’

’ den[t] = den[t] + numM[n,t];’

’ }’

’ for (n in 1:N){’

’ M[n,t] = numM[n,t]/ den[t];’

’ }’

’}’

’ F = UTree*M; ‘

’}’

’ model {’

’ for(t in 1:T) {’

’ for (n in 1:N){’

’ mu[n] ~ uniform(0, 0.3);’

’ delta_tau[n] ~ exponential(2);’

’ Fhat[n,t] ~ normal(F[n,t],0.01);’

’ }’

’ }’

’}’

};

%%

Fdata =

[1.0000,1.0000,1.0000,1.0000,1.0000,1.0000,1.0000,1.0000,1.0000,1.0000,1.0000

; 0, 0, 0, 0, 0, 0, 0, 0, 0,1.0000, 0

; 0,

0,0.0500,0.2000,0.2600,0.1000,0.2600,0.6000,0.6600,0.8800,0.9700

; 0,0.2800, 0, 0,0.7100,0.9100,0.7600,0.4600,0.4000,0.1300,0.0400

; 0, 0, 0, 0, 0,0.2500,0.3900,0.2200,0.4000,0.1400,0.0300

; 0, 0, 0, 0, 0, 0,0.0800,0.1200,0.1200,0.3400,0.400];

%%

Treedata = [ 0 0 0 0 0 0

;0 0 0 1 1 1

;0 0 0 0 0 0

;1 0 0 0 0 0

;0 0 0 0 0 0

;0 0 1 0 0 0];

UTreedata = [ 1 0 0 0 0 0

;1 1 1 1 1 1

; 0 0 1 0 0 0

; 1 0 0 1 0 0

; 0 0 0 0 1 0

; 0 0 1 0 0 1];

%%

T = size(Fdata,2);

N = size(Fdata,1);

%% initial values

delta_tau_init = ones(N,1);

mu_init = zeros(N,1);

evol_tree_dat = struct(‘T’,T,‘N’,N,‘Fhat’,Fdata,‘UTree’,UTreedata);

evol_tree_initial_values = struct(‘delta_tau’,delta_tau_init,‘mu’,mu_init);

%%

model =

StanModel(‘init’,evol_tree_initial_values,‘file_overwrite’,true,‘verbose’,true,‘method’,‘sample’,‘model_code’,evol_tree_code,‘data’,evol_tree_dat);

%%

model.compile();

fit_vb = model.vb();

print(fit_vb);

Thank you.