Next, we implemented mash modeling using mashr (v.0.2.57)21 for each feature using the limma-voom-extracted effect sizes and s.e. across brain regions. We learned the correlation structure across the brain regions and used all features as an unbiased representation of the results to account for overlapping samples. After this, we calculated the canonical covariances. A strong set of features was determined condition by condition using mash_1by1; data-driven covariance was calculated with the strong set of features. Once calculated, we fitted the mash model to the full set of features and computed the posterior summaries for all features. Features were considered significant if they had an LFSR < 0.05.