Puzzled about the array ordering of LKJCC corr as xarray object

Hey Oriol, thanks for coming back to it - much appreciated!

I don’t think that was a wild goose chase though - I still don’t understand why the PearsonR on b0 is flip-rotated from the correlations in lkjcc_corr

In any case, a simplified, pooled version looks a bit like the following… I hope this helps (and I hope moreover that it’s not fatally flawed since it comes from a couple of papers and seems to work for me so far! Haha).

This is a parametric approximation of a claims development curve where the cumulative sum lr_csum_ci of claims increases (generally speaking) over time t for a cohort of insurance policies.

We observe only the value of lr_csum_ci (here as a Lognormal), but decompose the median into two components:

  1. A development curve here implemented as a Fisk (log-logistic) CDF that transforms t: _fisk_cdf(gc_omega, gc_theta, t), ranges [0, 1] has support on t in (0, +inf)
  2. An expected loss ratio (ELR), has support on (0, +inf)

This is a useful decomposition for parameter interpretation / decision-making / setting priors according to business plan etc

So b0 is not directly observed - actually maybe I ought to non-center that, per the docs. Also note I don’t actually use corr_ (lkjcc_corr) for anything.

# handle data
y = pm.Data('y', obs['lr_csum_ci'].values, dims='obs_id')
t = pm.Data('t', obs['duration'].values + 1, dims='obs_id')

# correlated hyperpriors for GC and ELR intercept
sd_dist = pm.InverseGamma.dist(alpha=11, beta=10, shape=3)
chol, corr_, stds_ = pm.LKJCholeskyCov('lkjcc', n=3, eta=2, 
                              sd_dist=sd_dist, compute_corr=True)
b0 = pm.MvNormal('b0', mu=0, chol=chol, dims='b0_names')

# growth curve (cumulative) log-logistic (Fisk) CDF
gc_omega = tt.exp(b0[0])  # omega (slope at median)
gc_theta = tt.exp(b0[1])  # theta (median)
# elr
elr = tt.exp(b0[2])

yhat_s = pm.InverseGamma('yhat_s', alpha=101, beta=10)
yhat_s_w = yhat_s / (1 + (t / 40))  # weighting toward later developments
_ = pm.Lognormal(
    'yhat',
    mu=tt.log(_fisk_cdf(gc_omega, gc_theta, t) * elr),
    sigma=yhat_s_w,
    observed=y,
    dims='obs_id',
)

Cheers!