123456789101112131415161718192021222324252627282930313233343536373839 |
- from pymer4.models import Lmer
- def hierarchical_model_regions(dataframe):
- # overall p-value
- df_reduced = dataframe[['logepsilon', 'patient', 'binning_focus', 'region']]
- model = Lmer('logepsilon ~ binning_focus*region + (1|patient)',
- data=df_reduced)
- model.fit()
- anova_result = model.anova(force_orthogonal=True)
- global_pvals = anova_result['P-val']
- # separate p-values
- pvals_regions = dict()
- for reg in dataframe.region.unique():
- print(reg)
- df_region = dataframe[dataframe['region'] == reg]
- df_reduced = df_region[['logepsilon', 'patient', 'binning_focus']]
- # for formulas, check
- # https://stats.stackexchange.com/questions/13166/rs-lmer-cheat-sheet
- model = Lmer('logepsilon ~ binning_focus + (1|patient)',
- data=df_reduced)
- model.fit()
- anova_result = model.anova(force_orthogonal=True)
- region_pval = anova_result['P-val'].item()
- pvals_regions[reg] = region_pval
- return global_pvals, pvals_regions
- def hierarchical_model_hemispheres(dataframe):
- df_reduced = dataframe[['logepsilon', 'patient', 'binning_focus']]
- model = Lmer('logepsilon ~ binning_focus + (1|patient)',
- data=df_reduced)
- model.fit()
- anova_result = model.anova(force_orthogonal=True)
- pval = anova_result['P-val'].item()
- return pval
|