I have measurement data of a real experiment and a simulation modelling the same experiment. I now want to fit the simulation’s parameters, such that the simulation accurately resembles the real experiment.
The measurement data is given in form of multiple noisy time series. For simplicity let’s say the measurements are 10 time series of 100 entries each. The simulation has two parameters P,Q. For each parameter combination (P=p,Q=q) the simulation returns one time series of 100 entries.
I am struggling to tell pymc3 that I want only one single distribution for each of the two parameters, while the function evaluations are of length 100 and I have multiple reference time-series. Can someone please give me a hint or a quick demo on how to set up this scenario?
I’m not sure if I understood correctly. Maybe more information or some code would help.
Nevertheless, from what I understood you have 10 time series with 100 entries each. The process that is generating each time series is dependent on two parameters. Trying to be less vague, let’s suppose that we are talking about an AR(1) process that is generating each time series (we have two parameters, theta and tau).
Then,
N, s = 100, 10 # number of points per series and number of series
x = np.zeros((N,s))
true_thetas = np.repeat(np.array([[0.7]]),s).reshape(1,s) + np.random.normal(0, 0.2, s)
true_taus = stats.halfnorm.rvs(0.1, size=s)
x[0]=0.0
for i in range(1, N):
x[i] = true_thetas[0] * (x[i-1]) + stats.norm.rvs(true_center, np.sqrt(1/true_tau), s)
with pm.Model() as m_ar:
theta = pm.Normal('theta', 0, 1, shape=s) # The shape here is stating that there is only 1 theta for each series
tau = pm.HalfNormal('tau', sd=1, shape=s) # The shape here is stating that there is only 1 tau for each series
# AR Process
η_t = pm.AR1("η_t", k=theta, tau_e=tau, observed=x) # Here we are using the observed data, so it will evaluate for the 100 parameters of each series
trace = pm.sample()
