I want to estimate the parameters l and nu using the stan optimizing for the following stan model:
data {
int<lower=1> Mn; // number of nodes
int<lower=1> Mm; // number of modes
int<lower=0> N; // number of samples
row_vector[Mm] y[N];
matrix[Mm,3Mn] Umodes; // each row is a mode and each column is a node coordinate
row_vector<lower=0>[Mm] singval;
matrix<lower=0>[Mn,Mn] dij;
}
parameters {
row_vector[Mm] lambdaHat;
real var1;
real var2;
}
transformed parameters {
cov_matrix[Mm] sigmareduced;
{
real l=exp(var1);
real nusq=exp(var2);
matrix[Mn,Mn] pcordsigma=exp(-(1.0/l)dij); // Full covariance for a single coordinate
matrix[Mm,Mn] V1=Umodes[1:Mm,1:Mn]pcordsigma;
matrix[Mm,Mn] V2=Umodes[1:Mm,(Mn+1):(2Mn)]pcordsigma;
matrix[Mm,Mn] V3=Umodes[1:Mm,(2Mn+1):(3Mn)]pcordsigma;
sigmareduced=nusq(V1((Umodes[1:Mm,1:Mn])‘)+V2*((Umodes[1:Mm,(Mn+1):(2Mn)])')+V3((Umodes[1:Mm,(2Mn+1):(3Mn)])’));
sigmareduced=0.5*(sigmareduced+sigmareduced’);
}
}
model {
lambdaHat ~ std_normal();
//Likelihood
y ~ multi_normal(lambdaHat .* singval, sigmareduced);
}
generated quantities{
real nu=exp(var2/2);
real l=exp(var1);
}
Now, to improve the process, I would like to pass an initial value for the parameters. I saw that pystan allows specifying an initialisation using init followed by an init function that returns a dictionary.
Hence, my questions are:
-
Should I pass an initial value for all the entries in the section parameters (i.e., var1,var2 and lambdaHat)?
-
Since the initialisation depends on a parameter (that I will call nmodes) that furnishes the size of lambdaHat, is it correct to define the function as follows:
def initfun(nmodes):
l = math.exp(np.random.normal(loc=math.log(108.6),scale=0.2 ) )
nu = math.exp(np.random.normal(loc=math.log(113.7),scale=0.04 ) )
return dict(lambdaHat=np.random.normal(loc=0.0,scale=0.04,size=(nmodes)),var1=math.log(l), var2=math.log(2.0*nu) )
and pass to pystan as follows:
optim = sm.optimizing(data=inputdata,init=initfun(nmodes) )
When I tried, pystan did not raise any error;
Is this the correct way to specify an initialization?
Here initfun has one input argument and the input parameter isn’t the chain index.
Thanks,
Cesare