SSIMmap manuscript reproduction (two-sided raster update, 28 October 2026)
========================================================

This directory is installed with the package and located by:
  system.file("reproduction", package = "SSIMmap")
All analytical input values and geometries are included in the package datasets:
Toronto_SSIM, fwi_0816_bc, fwi_0818_bc and fwi_1111_bc. No downloads or absolute
computer-specific paths are needed. fwi_1101_bc is a legacy November 1 dataset;
it is deliberately NOT used for the revised manuscript's November 11 analysis.

Run in R >= 3.6.0 after installing SSIMmap >= 0.4.3 containing these files:

  source(system.file("reproduction", "reproduce_manuscript.R",
                     package = "SSIMmap"))
  reproduce_manuscript("ssim-observed", full = FALSE)
  reproduce_manuscript("ssim-full", full = TRUE)

Use a new writable output directory each time. All generated results are written
there, never in the installed package. No library is installed or modified by the
script. Public entry points restore the caller's RNG kind and random state on exit. The output VERIFICATION.txt explains the completed checks and limits.
Fast mode takes seconds to minutes. Full mode includes 6,000 FWI permutations;
the 3x3 local tests now replay those draws in a second pass to calculate
two-sided tail counts without retaining all null rasters. Runtime and memory use
depend on hardware and spatial-library versions.

Full mode generates:
- Table1_reproduced.csv: 24 means and all global two-sided p and BH q values;
- Table2_reproduced.csv under Toronto_tests: 16 means and two-sided global p;
- two 3 x 3 FWI local rasters including observed, p, BH q and significance layers;
- all four Toronto local results with observed values and fixed-seed p/q for SSIM, SIM, SIV and SIP;
- Figure 6 bandwidths 87 and 91 and all bias/variance curve values;
- numerical data for Figures 3-7, session information and comparison reports.

The two tables use raw (unadjusted) global p values for significance stars.
The raster local tests apply BH separately to each component's valid local cells.
The polygon local tests apply BH separately to each statistic's 279 local p values
for each comparison. Both data types use q <= alpha (or p <= alpha without FDR).
The polygon Figure 7 significance threshold is q < 0.01. The general sig column
returned by the package uses alpha = 0.05 and is NOT the Figure 7 classification.

Scientific settings and reference values
---------------------------------------
FWI: raw input values (transform="none"), k1=.01, k2=.03, map A 16 August 2023;
map B 18 August or 11 November 2023; window sides 3,5,7 (half-windows w=1,2,3).
R=1000; map A fixed, map B shuffled over common nonmissing input cells; two-sided
p=(1+sum(abs(null-mean(null)) >= abs(observed-mean(null))))/1001.
The center is specific to each global component or local cell/component. Global q adjusts 4 component p
values within each comparison/window. Seeds: summer_summer 202610271/272/273;
summer_fall 202610281/282/283 respectively.
Toronto: 279 ADA polygons; percentile transform; k1=.01,k2=.03; Gaussian adaptive
k-nearest-neighbor weights; bandwidths 17/87 for CIMD-PP and 17/91 for CIMD-Commute.
Both maps are independently permuted, R=1000; two-sided tests are centered on the
permutation mean. Seeds, respectively: 2026102701/2702/2703/2704. Figure 7 uses the
same seed as the corresponding Table 2 midpoint comparison. RNG kind is
Mersenne-Twister/Inversion/Rejection. Fixed-seed Figure 7 q<.01 counts are 279/124.

The expected/ CSV files contain the current fixed-seed references. Table1.csv and
Figure4_local_inference.csv use the 0.4.2 two-sided raster tests. The former
upper-tail references are retained under historical_0.4.1/ for comparison only.
analysis_manifest.json separates the current settings from the historical
0.4.0 audit-archive provenance; those external paths are not installed files.
Comparison tolerances are 1e-10 for table/curve numeric values and 1e-12 for local
polygon values. These allow floating-point noise but are much tighter than the
four decimals printed in the manuscript. Unexpected differences stop the script.

The original Figure 7 commuting local q values came from an unrecorded random seed.
They had 118 q<.01 areas; the new fixed-seed analysis has 124, with six changed
classifications. The revised analysis uses the fixed-seed results. The
archived legacy values remain only for a transparent before/after audit; they are
never substituted for freshly calculated p/q values.

Implementation verification
---------------------------
All three functions default to transform="none" in version 0.4.2. The reproduction
continues to specify none for FWI and percentile for Toronto, so all observed
values, polygon inference, and bandwidth results are unchanged. Raster p/q values
are recalculated using the permutation-mean-centered two-sided definition.
For 3x3 FWI, a private function copy records global null means during the first
pass of the exported local two-sided algorithm. Before the full run, the script
verifies exact equality of local outputs, global null means, and final RNG state
against the exports for 10 seeded draws. Other FWI cases use the exported global
function directly. Polygon verification remains unchanged. No installed function
or package file is modified by the reproduction script.

Optional independent FWI jobs / resume
-------------------------------------
Each of the six cases can run in a separate R session, using the same output root:
  source(system.file("reproduction", "reproduce_manuscript.R", package="SSIMmap"))
  run_fwi_case("ssim-full", "summer_summer", 1L)
Repeat with w=2L,3L and scenario="summer_fall". Do not run the same case twice in
one output root. After ALL six jobs have completed, assemble and verify them:
  reproduce_manuscript("ssim-full", full=TRUE, reuse_completed_fwi=TRUE)
Do not run the fast mode in that same output root before this assembly step.
Completed case reuse is explicit; every reused summary is rechecked against the
expected table. Outputs should be retained together with their sessionInfo files.

Scope and limits
----------------
This reproduces the analysis from the packaged analytical inputs, including the
current known historical Toronto aggregation. It does not derive Toronto_SSIM
anew from every underlying external Census/CIMD/MSDI source record.
Commute is the unweighted mean of stored DA percentages over all DA geometries
intersecting an ADA, including DA geometries that touch the ADA boundary. Each DA
percentage uses the number of people commuting to a different CSD AND CD within
the same province/territory, divided by TOTAL DA POPULATION, then rounded to two
decimals. It is not the ADA's pooled proportion, not all outside-CSD commuters,
not a household-population denominator and not restricted to member DAs. Sixteen
publicly suppressed source commuting cells are stored as 0 in the historical input;
their undisclosed counts cannot be recovered from the public Census. Existing
values are preserved and these limitations are stated in the dataset documentation.

The exports under figure_data reproduce scientific values and geometry. Original
publication TIFF typography, QGIS layout, fonts, basemaps, and CMYK profiles are a
separate rendering step, not a numerical claim of this script. Where the optional
figures/render_figures.R helper is present, see its documentation for rendering.
