# Preregistered analysis plan: Small-study inflation of whole-brain effect sizes in published unthresholded group fMRI maps (NeuroVault)

Operator: Quiet Replication. Preregistered before any outcome metric was computed on the frozen set (see "Pilot disclosure").

## Question and claim under test
Across published, group-level, unthresholded fMRI T and Z statistic maps shared on NeuroVault, do maps from smaller-N studies show larger **whole-brain** standardized effect sizes than maps from larger-N studies?
Prior evidence for small-study inflation in fMRI comes from reported peaks or thresholded clusters, which are inflated by selection (winner's curse) regardless of publication or analytic bias. A whole-brain, unthresholded metric is not subject to peak selection, so a negative slope on log N here cannot be explained by peak selection alone. A clean null is also informative and will be reported as a negative result.

## Sampling frame (frozen before any map was analysed)
- Source: NeuroVault public REST API (https://neurovault.org/api/), crawled 2026-10-08 16:12–16:50 PT with a curl user agent. All 17,948 public collections were paged; for every collection with a non-empty publication `DOI` and ≥1 image, all of its images were paged (918 collections, 173,015 image records). Raw API pages are kept (gzipped) for audit.
- Inclusion filter (applied in this order, code `code/freeze.py`):
  1. collection has a publication DOI (field `DOI`; preprint-only DOIs not counted);
  2. `image_type` = statistic_map; 3. `modality` = fMRI-BOLD; 4. `analysis_level` = group;
  5. `map_type` ∈ {T map, Z map}; 6. `is_thresholded` = False; 7. `not_mni` = False;
  8. volumetric (`data_origin` volume or empty); 9. species homo sapiens (or empty); 10. file is .nii/.nii.gz;
  11–12. `number_of_subjects` present, integer, N ≥ 10;
  13. collection does not declare `group_comparison` = True;
  14. image text (name, description, contrast definition, Cognitive Atlas contrast, file name; `_` and `.` read as spaces) does not match the between-group regex `GROUP_RE`;
  15. image text does not match the correlation/covariate regex `CORR_RE`;
  16. the map is a task-activation contrast: Cognitive Atlas paradigm is not a resting-state paradigm, and image text does not match `NONTASK_RE` (ICA/components, seed/PPI/functional connectivity, resting state, ALFF/ReHo/degree centrality, F-tests, MVPA/decoding/searchlight/RSA, single-subject names);
  17. for T maps whose metadata gives degrees of freedom: exclude if df > N−1+0.5 or df < 0.5·(N−1) (N inconsistent with df).
  All regexes are verbatim in `code/freeze.py`. Text classification is "as far as the metadata says": misclassified maps are expected to remain at a low rate and are not hand-edited.
- Frozen list: `included_maps.csv` (1,750 maps from 258 collections), **SHA-256 `9c2b3b7e37a2cffea63b9684cc2a7cc55ca088976a22ad6cfb6c958d3acdad86`**. Columns: image_id, collection_id, file_url, N, map_type, df_meta, fwhm_mm (image `smoothness_fwhm`, else collection `smoothing_fwhm`, else empty), upload_year (image `add_date` year), task (Cognitive Atlas paradigm name, else `unknown`), task_id, doi, name, contrast_definition. Exclusion counts at each step: `freeze_summary.json`.
- N distribution of the frozen list: min 10, 10th pct 19, quartiles 23 / 32 / 77, 90th pct 244, max 1,369; maps with N 10–20: 298, 21–30: 525, 31–50: 349, 51–100: 235, 101–200: 143, 201–500: 194, >500: 6. 1,374 T maps, 376 Z maps; median 4 maps per collection (max 176; 10 collections contribute >20 maps); collection-median N quartiles 21 / 29 / 49. Paradigm `unknown` for 16% of maps; FWHM missing for 74%; upload years 2016–2025. Only one paradigm (go/no-go task, 11 collections) currently meets the S15 ≥10-collection threshold..

## Data acquisition (phase 2, after this preregistration is accepted)
Each listed map is downloaded once from its NeuroVault URL (`code/fetch_maps.py`); its SHA-256 is pinned at first fetch in `data/map_pins.csv`; maps that cannot be downloaded (HTTP error after retries) are excluded as `X1_download_failed` and counted. The reporting bundle vendors the map files so the Docker analysis runs offline; `code/analyze.py` refuses to run unless the frozen list and mask hashes match this plan, and checks every map against its pin.

## Mask and resampling (locked)
- Fixed 2 mm MNI152 (ICBM 2009c, nilearn-packaged) grey-matter mask: nilearn 0.11.1 `load_mni152_gm_mask(resolution=2, threshold=0.2, n_iter=2)`, 204,492 voxels, file SHA-256 `8e9037ed23d20fc495d63905ee75faf0c10203e62dd5a30b86a0b3edac409412` (`code/make_mask.py`).
- Each map is sampled at every mask-voxel centre by **nearest-neighbour** lookup through the image affine (no interpolation: linear/spline interpolation would shrink the voxelwise variance of t and bias the metric downward for coarser maps). 4D maps with one volume are squeezed; other 4D maps are excluded (X3).
- Valid voxels: finite and non-zero (zeros are treated as outside the analysed field).

## Outcome metric (locked)
For each map with N subjects and valid in-mask values v (t or z):
- T maps, ν = N − 1: δ̂² = ( mean(t²) · (ν−2)/ν − 1 ) / N.
- Z maps: δ̂² = ( mean(z²) − 1 ) / N.
Rationale: for a one-sample/paired t with noncentrality δ√N, E[t²] = ν(1+Nδ²)/(ν−2); for z, E[z²] = 1 + Nδ². So δ̂² is a moment estimator of the mean squared standardized (Cohen's d-scale) effect over the grey-matter mask; spatial correlation does not bias a mean.
- δ̂_rms = sqrt(max(δ̂², 1e-4)); outcome y = log δ̂_rms.
- **Non-positive values:** maps with δ̂² ≤ 1e-4 (including δ̂² ≤ 0) are kept and set to the floor (δ̂_rms = 0.01, y = log 0.01). Their count is reported. Secondary S1 drops δ̂² ≤ 0 maps instead; S2 uses floor 1e-3.

## Quality control (locked; applied after download, before modelling)
- X3: unreadable file, not 3D, or degenerate affine.
- X4 coverage: fewer than 70% of mask voxels valid (or < 1000 valid voxels).
- **X5 core test ("least-extreme voxels inconsistent with a t/z distribution with the stated N"):** let q10 = 10th percentile of |v| over valid voxels and c10 = the null 10th percentile of |T|, i.e. t_{N−1}^{-1}(0.55) for T maps or Φ^{-1}(0.55) = 0.125661 for Z maps. Exclude if **q10 < 0.5 · c10**. Justification: the folded noncentral t (or normal) is stochastically increasing in |noncentrality|, so for a genuine t/z map every voxel's |v| is stochastically at least as large as under the null, and the empirical 10th percentile cannot fall far below c10 regardless of the true effects; the factor 0.5 allows for spatial correlation. This catches maps that are not standardized statistics (betas/contrast estimates, percent signal change, -log10 p, F-like or post-hoc rescaled maps) labelled as T/Z.
- **X6 extreme-tail plausibility:** exclude if the 99.9th percentile of |v| > 4·√N + 8 (would imply voxel effects of d ≳ 4, which indicates a mislabelled map or an understated N).
- X7 content duplicates: among QC-passing maps, identical masked content (values rounded to 3 dp + validity pattern) → keep the lowest image_id (removes the same group map uploaded to several collections).
- Honest limitation (design refinement, see below): the core test detects wrong map *type*; it cannot detect a wrong N, because the low quantiles of t_{ν} barely depend on ν. Wrong N is screened by the metadata df check at freeze (step 16) and the X6 tail check, and the S14 sensitivity reports results without X5/X6.

## Primary model and statistic (locked)
Linear mixed model (statsmodels 0.14.4 `MixedLM`, REML), one row per analysed map:
y = log δ̂_rms ~ β·log N + map type (T/Z) + FWHM (mm; 0 when missing) + FWHM-missing indicator + upload-year fixed effects + paradigm fixed effects (Cognitive Atlas task, else `unknown`) + random intercept per collection.
Optimizers tried in order lbfgs, bfgs, powell, nm (maxiter 5000); the first converged fit is used. Linearly dependent design columns are dropped deterministically (earliest kept; log N always retained). 95% CI: Wald normal.

## Decision rule (locked)
- **Support** (smaller-N maps show larger whole-brain effects): β ≤ −0.10 **and** 95% CI upper bound < 0.
- **Refute** (no meaningful small-study inflation): 95% CI lower bound > −0.10.
- **Inconclusive** otherwise.
(β = −0.10 means a 10-fold smaller N goes with ~26% larger rms effect.)

## Preregistered secondary analyses (descriptive; not part of the decision rule)
S1 drop δ̂² ≤ 0 instead of flooring; S2 floor 1e-3; S3 one map per collection (seeded random pick, seed 20261008; OLS, HC3) and S3b the median and 2.5/97.5 percentiles of β over 200 such draws; S4 exclude collections contributing > 20 maps, and > 5 maps; S5 collection-level OLS of mean y on log(median N) + fraction T maps; S6 T maps only; S7 Z maps only; S8 N ≤ 100; S9 add a contrast-vs-baseline indicator (regex `baseline|\brest(ing)?\b|fixation|\bfix\b|implicit|\bnull\b|\boff\b` on name + contrast definition) because large-N consortium uploads may favour condition-difference contrasts while small studies upload task-vs-baseline maps; S10 OLS with collection-clustered SEs; S11 no paradigm FE; S12 known paradigms only; S13 keep content duplicates; S14 QC X3/X4 only (keep X5/X6 failures).
**S15 leave-collection-out ROI Cohen's d:** for paradigms with ≥ 10 distinct collections among analysed maps, for each map build the consensus d-map (equal-weight mean over the *other* collections of their collection-mean v/√N maps), take the top 1% of mask voxels by |consensus| (2,045 voxels; stable sort), and compute d_ROI = mean over valid ROI voxels of sign(consensus)·v / √N (signed, so contrast-direction mismatches add noise rather than bias). Model d_ROI with the primary covariates and collection random intercept; report β on log N. Because the ROI is defined without the map's own collection, it carries no selection bias.
Each secondary also reports whether it would meet the primary support/refute thresholds (`direction_rule`).

## Interpretation limits (stated in advance)
A negative β is consistent with publication/analytic bias but can also arise from genuine heterogeneity (small studies choosing stronger paradigms or contrasts, more homogeneous samples, different pipelines); paradigm FE, S9 and S15 address this only partly. The claim is about published maps shared on NeuroVault, not all fMRI studies.

## Pilot disclosure
Before this preregistration, the following map files were downloaded and processed only to test code: five listed maps drawn with numpy default_rng(20261008) from the frozen list — image ids 24834, 42337, 43172, 50031, 869463 (per-map pilot outputs in `pilot/per_map_pilot.csv`; one failed X4 coverage, four produced δ̂ values; no model was or can be fit on them). Separately, during scouting, image 601944 (collection 11654, not in the frozen list: no publication DOI) was downloaded to inspect the file format. The scout phase also downloaded image 32015 (single-subject, not eligible) to inspect file formats. No outcome metric or model was computed on any other map. The analysis code was otherwise tested end-to-end only on synthetic simulated maps (`code/tests/make_synthetic.py`; recovers the planted slope).

## Reproducibility
Python 3.12 (public.ecr.aws/docker/library/python:3.12-slim), pinned deps in `env/requirements.txt`; `code/run` executes `code/analyze.py` offline with single-threaded BLAS; all reported floats rounded to 6 significant digits / 6 decimals; seed 20261008.

## Design changes relative to the scouting proposal
1. Crawl: instead of paging all ~715k image records, all 17,948 collections were paged and then every image of every DOI-bearing collection (the DOI is a collection-level field, so the included set is identical); completeness verified (173,015 records = sum of per-collection counts).
2. Added inclusion step 16 (task-activation maps only). Inspection of metadata (no map values) showed ~900 resting-state, ICA, seed/PPI connectivity, ALFF/ReHo, F-test and decoding maps. Their T/Z values are not Cohen's-d-scale activation contrasts (ICA z-scores are not test statistics at all), so the δ̂ formula does not apply. Correlation and between-group regexes were also widened (e.g. `corr`, `relation with`, `group vs group`, `BMI`) after reading map names. All of this was decided before any outcome was computed beyond the 5-map pilot.
3. QC core test is one-sided (low quantile) and justified by stochastic ordering of the folded noncentral t; it detects wrong map *type* but cannot detect wrong N, so a metadata df/N consistency check (step 17) and an extreme-tail plausibility check (X6) were added. Understated N is the error that would spuriously favour the hypothesis.
4. Nearest-neighbour resampling, not interpolation (interpolation shrinks t variance and would bias δ̂ by voxel size).
5. Content de-duplication (X7) so the same group map uploaded to several collections is counted once.
6. Added S9, a contrast-vs-baseline covariate, as a confound check, and stated the heterogeneity caveat in advance.
7. The ROI-d secondary uses a signed, leave-collection-out consensus ROI (top 1%), so neither selection nor folding bias favours the hypothesis.

## Report-by
2026-10-22