Rhat Nan Metropolis


#1

I’m trying to see whether there was a switch point in Syrian conflict casualties (2011-2014, grouped by casualties per week).

with pm.Model() as model:
    alpha = 1.0/a.casualties.mean() 
    lambda_1 = pm.Exponential("lambda_1", alpha)    
    lambda_2 = pm.Exponential("lambda_2", alpha)    
    tau = pm.DiscreteUniform("tau", lower=0, upper=count_conflict_time-1)#switchpoint
   
    idx = np.arange(count_conflict_time) # Index
    lambda_ = pm.math.switch(tau > idx, lambda_1, lambda_2)
    observation = pm.Poisson("obs", lambda_, observed=a.casualties)

with model:
    trace = pm.sample(10000, tune=5000)

pymc3 automatically chooses the right sampler:

CompoundStep
NUTS: [lambda_2_log__, lambda_1_log__]
Metropolis: [tau]

But when look at summary (pm.summary(trace)) I see Rhat is NaN:

     mean    sd 	mc_error 	hpd_2.5 	hpd_97.5 	n_eff 	Rhat
 tau  78.0, 0.0      0.0,      78.0,         78.0,      20000.0, NaN

And by plotting trace we see:

No matter what I do, I always get scale reduction factors (rhat) = NaN. Does this mean the tau didn’t converge? How can this be corrected?


#2

As you can see from the trace, tau is completely stuck. This is likely due to the inefficient of random walk Metropolis in high dimension: all proposals are rejected.

While Discrete switch point model could be sampled in PyMC3 after some careful tuning, our recommendation is to rewrite them as a continuous model. You can have a look of this notebook to see how to do it:


#3

thank you, junpenglao. When I try to run your code above, I get:

Can’t get attribute ‘Composed’ on <module ‘main’ (built-in)>

I’ve been trying for hours to fix that to no avail. Do you have any idea what is wrong?

This suggestion says to write the function into a different file, but it didn’t seem to work.