`extract_components_from_idata` doesn't handle multivariate time series

[AlexAndorra]

There is a bug in extract_components_from_idata when working with multivariate time series. The following model will show it:

ll = pmss.structural.LevelTrendComponent(
    name="trend", order=2, innovations_order=[0, 1], observed_state_names=presidents # list of 5 names for 5 multivar TS
)
ar = pmss.structural.AutoregressiveComponent(name="ar", order=1, observed_state_names=presidents)
noise = pmss.structural.MeasurementError(name="noise", observed_state_names=presidents)
structural_mod = (ll + ar + noise).build()

Once fitted, post_pred = structural_mod.sample_conditional_posterior(idata) returns 15 hidden states, as it should:

array(['level[chirac2]', 'trend[chirac2]', 'level[sarkozy]', 'trend[sarkozy]',
       'level[hollande]', 'trend[hollande]', 'level[macron]', 'trend[macron]',
       'level[macron2]', 'trend[macron2]', 'L1[chirac2]', 'L1[sarkozy]',
       'L1[hollande]', 'L1[macron]', 'L1[macron2]'], dtype='<U15')

but then, structural_mod.extract_components_from_idata(post_pred) only returns one AR instead of 5:

array(['trend[level[chirac2]]', 'trend[trend[chirac2]]',
       'trend[level[sarkozy]]', 'trend[trend[sarkozy]]',
       'trend[level[hollande]]', 'trend[trend[hollande]]',
       'trend[level[macron]]', 'trend[trend[macron]]', 'trend[level[macron2]]',
       'trend[trend[macron2]]', 'ar'], dtype='<U22')

[AlexAndorra]

I think I understand what's happening: Autoregressive, Cycle and Seasonality always use combine_hidden_states=True, so _hidden_states_from_data sums all states into a single output state, which is incorrect for multivariate models.

For this model, changing the logic to combine_hidden_states=k_endog == 1 in all those components' super.__init__ call works, returning:

array(['trend[level[chirac2]]', 'trend[trend[chirac2]]',
       'trend[level[sarkozy]]', 'trend[trend[sarkozy]]',
       'trend[level[hollande]]', 'trend[trend[hollande]]',
       'trend[level[macron]]', 'trend[trend[macron]]', 'trend[level[macron2]]',
       'trend[trend[macron2]]', 'ar[L1[chirac2]]', 'ar[L1[sarkozy]]',
       'ar[L1[hollande]]', 'ar[L1[macron]]', 'ar[L1[macron2]]']

from structural_mod.extract_components_from_idata(post_pred).

But, IIUC, I need to find something more refined, so that if for instance we have an AR2 instead of AR1, we get back:

array(['trend[level[chirac2]]', 'trend[trend[chirac2]]',
       'trend[level[sarkozy]]', 'trend[trend[sarkozy]]',
       'trend[level[hollande]]', 'trend[trend[hollande]]',
       'trend[level[macron]]', 'trend[trend[macron]]', 'trend[level[macron2]]',
       'trend[trend[macron2]]', 'ar[chirac2]', 'ar[sarkozy]',
       'ar[hollande]', 'ar[macron]', 'ar[macron2]']

My current fix would not work and would just return ..., 'ar[L1[chirac2]]', 'ar[L2[chirac2]]', 'ar[L1[sarkozy]]', 'ar[L2[sarkozy]]', ..., basically duplicating what sample_conditional_posterior is doing (cc @jessegrabowski , in case I missed something)

[jessegrabowski]

When combine_hidden_states is true, the summation is done using the obs_state_idx provided by each component. See here.

Only one AR component should be returned, because only one state is relevant. For AR(3), the component is:

$$ \begin{bmatrix}x_t \\ x_{t-1} \\ x_{t-2} \end{bmatrix} = \begin{bmatrix} \rho_1 & \rho_2 & \rho_3 \\ 1 & 0 & 0 \\ 0 & 1 & 0 \end{bmatrix} \begin{bmatrix}x_{t-1} \\ x_{t-2} \\ x_{t-3} \end{bmatrix} $$

As you can see, if you look at the 1st state, you have all the information, the other states are just carrying along bookkeeping information.

[jessegrabowski]

I don't believe there is a bug here. There seems to be some confusion about what extract_components_from_idata does. It does not return all of the hidden states. That is the job of sample_xxx_posterior. The purpose of extract_components_from_idata is to perform time series decomposition. It should return data = level + seasonal + cycle + ar.

This is not the same as simply summing all the states. If you look at the AR example above, you will see summing everything results in $\rho_1 x_{t-1} + \rho_2 x_{t-2} + \rho_3 x_{t-3} + x_{t-2} + x_{t-3}$, which does not correspond to anything in the model.

[AlexAndorra]

But shouldn’t it return one AR component *per time series*? If my model has 5 AR components (one for each observed states), I’m expecting back 5 AR components from extract_components_from_idata , not just one as is currently the case — I don’t see why the observed states should be squeezed into one. Or am I missing something? Il giorno sab 19 lug 2025 alle 21:33 Jesse Grabowski < ***@***.***> ha scritto:

…

*jessegrabowski* left a comment (pymc-devs/pymc-extras#544) <#544 (comment)> I don't believe there is a bug here. There seems to be some confusion about what extract_components_from_idata does. It does *not* return all of the hidden states. That is the job of sample_xxx_posterior. The purpose of extract_components_from_idata is to perform time series decomposition <https://otexts.com/fpp3/components.html>. It should return data = level + seasonal + cycle + ar. This is not the same as simply summing all the states. If you look at the AR example above, you will see summing everything results in $\rho_1 x_{t-1} + \rho_2 x_{t-2} + \rho_3 x_{t-3} + x_{t-2} + x_{t-3}$, which does not correspond to anything in the model. — Reply to this email directly, view it on GitHub <#544 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/AHIJMTC2YZ53YTIXXMC2HH33JLWYBAVCNFSM6AAAAACB5ALXB2VHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZTAOJSGY4DCNRSGI> . You are receiving this because you were assigned.Message ID: ***@***.***>

[jessegrabowski]

Ah I missed the plot. Yes, this is a bug.