Error running linear mixed effects in the FMT tutorial

Hello,

I've started using (and learning) QIIME2 a couple of weeks ago and I'm still trying to run the tutorials on my own. I run most of the FMT tutorial by now, but I'm trying to answer some of the questions and I'm having issues with the q2-longitudinal plugin and linear mixed effects command.

I've tried the following code:

qiime longitudinal linear-mixed-effects \
  --m-metadata-file filtered-sample-metadata.tsv \
  --m-metadata-file core-metrics-results/faith_pd_vector.qza \
  --p-metric faith_pd \
  --p-group-columns sample-type,treatment-group,administration-route,gender \
  --p-state-column week \
  --p-random-effects age,weight,height,bmi \
  --p-individual-id-column subject-id \
  --o-visualization linear-mixed-effects-faith.qzv

And this error followed:

Plugin error from longitudinal:

  Shape mismatch between endog/exog and extra 2d arrays given to model.

Debug info has been saved to /tmp/qiime2-q2cli-err-e0qwm2q7.log

I also tried to remove the random effects parameters:

qiime longitudinal linear-mixed-effects \
  --m-metadata-file filtered-sample-metadata.tsv \
  --m-metadata-file core-metrics-results/faith_pd_vector.qza \
  --p-metric faith_pd \
  --p-group-columns sample-type,treatment-group,administration-route,gender \
  --p-state-column week \
  --p-individual-id-column subject-id \
  --o-visualization linear-mixed-effects-faith.qzv

But then there was a different error:

Plugin error from longitudinal:

  index 71 is out of bounds for axis 0 with size 71

Debug info has been saved to /tmp/qiime2-q2cli-err-whn65511.log

I've looked in the forum for any topic similar to this but I couldn't find. I would welcome any advice on how to better run these commands as everything is still quite challenging.

Thank you,

Hi @mgmartins
Can you run the verbose flag on one of those commands so we can get a better look at the issue?
:turtle:

Hi @cherman2,

I've just run the verbose flag like this:

qiime longitudinal linear-mixed-effects \
  --m-metadata-file filtered-sample-metadata.tsv \
  --m-metadata-file core-metrics-results/faith_pd_vector.qza \
  --p-metric faith_pd \
  --p-group-columns sample-type,treatment-group,administration-route,gender \
  --p-state-column week \
  --p-random-effects age,weight,height,bmi \
  --p-individual-id-column subject-id \
  --o-visualization linear-mixed-effects-faith.qzv \
  --verbose

And this was the output, which to be honest I don't understand

Traceback (most recent call last):
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/q2cli/commands.py", line 468, in __call__
    results = action(**arguments)
  File "<decorator-gen-521>", line 2, in linear_mixed_effects
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/qiime2/sdk/action.py", line 274, in bound_callable
    outputs = self._callable_executor_(
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/qiime2/sdk/action.py", line 558, in _callable_executor_
    ret_val = self._callable(output_dir=temp_dir, **view_args)
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/q2_longitudinal/_longitudinal.py", line 200, in linear_mixed_effects
    model_summary, model_results, model_fit, formula = _linear_effects(
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/q2_longitudinal/_utilities.py", line 353, in _linear_effects
    mlm = mixedlm(
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/statsmodels/regression/mixed_linear_model.py", line 1046, in from_formula
    mod = super(MixedLM, cls).from_formula(
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/statsmodels/base/model.py", line 229, in from_formula
    mod = cls(endog, exog, *args, **kwargs)
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/statsmodels/regression/mixed_linear_model.py", line 741, in __init__
    super(MixedLM, self).__init__(endog, exog, groups=groups,
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/statsmodels/base/model.py", line 270, in __init__
    super().__init__(endog, exog, **kwargs)
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/statsmodels/base/model.py", line 95, in __init__
    self.data = self._handle_data(endog, exog, missing, hasconst,
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/statsmodels/base/model.py", line 135, in _handle_data
    data = handle_data(endog, exog, missing, hasconst, **kwargs)
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/statsmodels/base/data.py", line 675, in handle_data
    return klass(endog, exog=exog, missing=missing, hasconst=hasconst,
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/statsmodels/base/data.py", line 72, in __init__
    arrays, nan_idx = self.handle_missing(endog, exog, missing,
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/statsmodels/base/data.py", line 257, in handle_missing
    raise ValueError("Shape mismatch between endog/exog "
ValueError: Shape mismatch between endog/exog and extra 2d arrays given to model.

Plugin error from longitudinal:

  Shape mismatch between endog/exog and extra 2d arrays given to model.

See above for debug info.

The output for the command without the random effects parameters was the following:

Traceback (most recent call last):
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/q2cli/commands.py", line 468, in __call__
    results = action(**arguments)
  File "<decorator-gen-521>", line 2, in linear_mixed_effects
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/qiime2/sdk/action.py", line 274, in bound_callable
    outputs = self._callable_executor_(
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/qiime2/sdk/action.py", line 558, in _callable_executor_
    ret_val = self._callable(output_dir=temp_dir, **view_args)
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/q2_longitudinal/_longitudinal.py", line 200, in linear_mixed_effects
    model_summary, model_results, model_fit, formula = _linear_effects(
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/q2_longitudinal/_utilities.py", line 353, in _linear_effects
    mlm = mixedlm(
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/statsmodels/regression/mixed_linear_model.py", line 1046, in from_formula
    mod = super(MixedLM, cls).from_formula(
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/statsmodels/base/model.py", line 229, in from_formula
    mod = cls(endog, exog, *args, **kwargs)
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/statsmodels/regression/mixed_linear_model.py", line 803, in __init__
    self.endog_li = self.group_list(self.endog)
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/statsmodels/regression/mixed_linear_model.py", line 1104, in group_list
    return [np.array(array[self.row_indices[k]])
  File "/home/marina/miniconda3/envs/qiime2-2023.5/lib/python3.8/site-packages/statsmodels/regression/mixed_linear_model.py", line 1104, in <listcomp>
    return [np.array(array[self.row_indices[k]])
IndexError: index 71 is out of bounds for axis 0 with size 71

Plugin error from longitudinal:

  index 71 is out of bounds for axis 0 with size 71

See above for debug info.

Hi @mgmartins,
I have seen this before but I can't quite remember what the issue was. I think it was something about sample ids in the metadata and in the diversity metric not completely lining up. (But I could be wrong)
Would you mind sending the metadata and the faith_pd_vector for me to look at?

Here are the metadata file I'm using, @cherman2, and I believe the faith_pd_vector you want, but I'm not sure. Let me know if it is another file

Thank you!

alpha-diversity.tsv (6.4 KB)

filtered-sample-metadata.tsv (10.7 KB)

@mgmartins,
It does seem like your IDs in those files your sent me don't entirely match up. I think your metadata has less sample ids then your faith_pd_vector you sent me. That could possibly be it?

Could you check to see if your IDs from your feature table and metadata are the same?

Also would you send me the qza of your faith_pd_vector.qza so that I can try to run this command on my end?

One more thing: can you link the tutorial you are referring to?

Sorry I don't have a direct answer but I will keep trying to help until we got these commands running!!

Hi @cherman2,

I don't know how to check the IDs in the feature table... One thing I did was to filter only treatment and control subjects to run the longitudinal analysis. Maybe it was not the best way, as I did it after running the q2-diversity core-metrics-phylogenetic but before doing the longitudinal analysis. Could this new table be the problem?

I'm trying to answer the Diversity analysis questions here, using as an example the linear mixed model of the longitudinal tutorial.

Thank you so much for your help, I'm really learning a lot!

faith_pd_vector.qza (98.3 KB)

Hi @mgmartins,
Thank you for your patience

Good thinking! This shouldn't be the problem because that new table is not an input to this command so this command has no idea that you filtered a table!

I am going to run these commands and see if I can figure this out! I will be back with an update.
:turtle:

Hi again!
So I was able to re-create your error.

I still think its because some sample ids are in the alpha diversity metric and not in the metadata.
Below is the list of samples that are missing from the metadata.

'8d2746d9.68e2.42e3.acfc.092603ad6ec1',
'6519888d.bb69.495b.b4b8.07d0d4333455',
'78ff0b74.3d7e.4ed5.8c2d.a5f77b793718',
'b986755d.b15a.4ad4.b6a2.0c7cdd9b8ced',
'2756472a.03f4.41da.8564.1d7538f949c6'

Hope that helps!
:turtle:

Hi @cherman2,

Thanks for all your help! Those five sample missing in the metadata are the donor samples I've filtered from the analysis, so maybe that is what is happening.

I've tried to run the same command just using the treatment-group in the --p-group-columns and it worked, so I do believe it is something I've done with the filtered metadata. Maybe if there are some samples missing, it is not able to group them by the other categorical variables.

qiime longitudinal linear-mixed-effects \
  --m-metadata-file filtered-sample-metadata.tsv \
  --m-metadata-file core-metrics-results/faith_pd_vector.qza \
  --p-metric faith_pd \
  --p-group-columns treatment-group \
  --p-state-column week \
  --p-random-effects gender,age,weight,height,bmi \
  --p-individual-id-column subject-id \
  --o-visualization linear-mixed-effects-faith.qzv
1 Like

This topic was automatically closed 31 days after the last reply. New replies are no longer allowed.