{"id":"2efca5f3-720a-4ffb-ab40-fef370e5c16a","arxiv_id":"2504.15586","paper_version":1,"verdict":"CONDITIONAL","confidence":"HIGH","novelty_score":6.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":3,"one_line_summary":"Joint leave-group-out scoring of spatial predictions improves Bayesian spatial model selection accuracy and reduces variability compared to pointwise scoring when spatial dependence is strong.","lead":"This paper tests whether scoring spatial cross-validation predictions as a joint block, rather than one point at a time, improves model selection for spatial data. It finds in simulations and two real datasets that joint scoring reduces variability and improves selection accuracy when spatial dependence is strong.","discovery_kind":"extension","skeptic_critique":{"model":"deepseek-v4-flash","headline":"The headline result rests on simulated CV scores computed with a Laplace approximation that is validated only for two kernel-model marginals, not for the SAR experiments or the held-out predictive densities; if the approximation biases joint and pointwise scores asymmetrically, the central claim…","rationale":"The paper's central claim is an empirical one: under strong spatial dependence, joint scoring yields more reliable model selection than pointwise scoring. The evidence for this is the simulation study, whose quantitative outputs (correctness, Z ratio) are computed from predictive densities obtained via the Appendix D Laplace approximation. The unusual element is that the approximation is applied to the joint posterior of parameters and held-out data, and the predictive density is evaluated at observed test values after optimizing over those test values. This is not a standard off-the-shelf Laplace approximation, and the validation provided is limited to two marginal posteriors from one data draw in the kernel experiments. The SAR experiments, which contain the headline rho threshold, are not validated at all. Because joint scoring uses the full predictive covariance matrix and pointwise scoring does not, any systematic error in the approximated predictive covariance—e.g., underestimation of variance under strong dependence—will affect the two scores differently. This makes the simulation result potentially an artifact of the approximation. The concern is concrete and testable: rerun a small slice of the SAR experiment with MCMC or exact Gaussian inference. If the conclusions survive, the paper's central claim is supported; if not, it is unsupported. The reader already issued a CONDITIONAL verdict, and this concern reinforces that condition rather than changing it. I do not see a more fundamental flaw: the experimental design, paired comparisons across test set sizes, and use of the Z ratio are reasonable, and the paper is honest about limitations such as the log score and the SAR parameterization.","tokens_in":26322,"tokens_out":9076,"duration_ms":94766,"concrete_test":"Run exact (or MCMC) inference for a focused subset of the Section 3.1.1 covariate selection design: rho*=0.99, ntest=1 and ntest=64, N=100 replications, using the same priors and candidate models. Compute correctness and Z ratio under pointwise and joint scoring and compare with Figure 2. If the joint advantage over pointwise disappears or reverses, the simulation result is an artifact of the Laplace approximation; if it persists, the central claim survives this check.","verdict_should_be":"UNCHANGED","load_bearing_attack":"The main quantitative support for the central claim is the simulation study of Section 3, and every simulated CV score is computed with the Laplace approximation described in Appendix D. The procedure in Eqs. (D.1)-(D.2) optimizes over the held-out test and buffer values before forming a Gaussian approximation from the Hessian. The only validation (Figures D.1, D.2) compares Laplace and MCMC marginal posteriors for a single data draw in the covariance-kernel experiments of Section 3.2. It does not validate the SAR experiments (Sections 3.1.1, 3.1.2) that establish the headline rho>0.75 result, nor does it validate the held-out predictive densities p(ytest|ytrain) that actually enter the joint and pointwise scores. Joint scores use the full predictive covariance, pointwise scores use only marginal variances, so an approximation error in the predictive covariance could change the two scoring methods asymmetrically. If the Laplace approximation becomes more biased as rho and block size grow, the apparent increase in Z ratio under joint scoring could be an artifact. The paper's own text acknowledges that Laplace approximation 'can be a crude posterior approximation,' which makes this the least secure link in the argument.","agreement_with_reader":"partial"},"referee_report":{"model":"deepseek-v4-flash","summary":"The paper studies block-wise cross-validation (CV) for Bayesian spatial models with Gaussian covariance structures, focusing on whether jointly evaluating the predictive performance of left-out blocks (via the joint log predictive density) is preferable to aggregating pointwise predictive scores. Through simulation experiments on a regular lattice, covering SAR-based covariate selection, network structure selection, and covariance kernel selection, the paper reports that joint scoring yields higher model selection accuracy and lower relative variability (measured by a Z ratio) than pointwise scoring, particularly when spatial dependence is strong (standardized rho > 0.75). Two real-data case studies (Australian child vaccination rates and Pennsylvania lung cancer incidence) illustrate the approach and produce results broadly consistent with the simulations. The paper acknowledges the limitations of using the Laplace approximation for computational feasibility and discusses the limited generalizability of the standardized-rho interpretation.","tokens_in":26580,"tokens_out":5602,"duration_ms":52579,"significance":"The paper addresses an under-explored aspect of spatial cross-validation: the scoring rule used to evaluate left-out blocks, as opposed to the more studied blocking design. If the empirical claim is reliable, the recommendation to use joint evaluation for spatially clustered CV would be useful for practitioners, as it extends findings from time-series settings (Cooper et al., 2024) to spatial models and to model components beyond covariates, such as network structure and covariance kernels. The simulation design is a strength: common data draws are used across test set sizes, multiple model selection tasks are considered, and the paper is transparent about the limitations of its simulation-based approach. However, the central claim currently rests on simulation results whose computational approximation is only weakly validated and whose statistical uncertainty is not quantified, so the strength of the conclusions exceeds what the evidence presently supports.","major_comments":[{"comment":"All simulation-based CV scores in Section 3 are computed with the Laplace approximation described in Appendix D, but the only validation against exact MCMC inference is for marginal posterior distributions of hyperparameters in the two covariance-kernel models (Section 3.2), for a single data draw. The SAR experiments (Sections 3.1.1 and 3.1.2) that support the headline rho > 0.75 claim are not validated, nor are the held-out predictive densities p(ytest_k | ytrain_k) that enter the joint and pointwise scores in Eqs. (2) and (3). This is a load-bearing gap because joint scores use the full predictive covariance, while pointwise scores use only marginal variances; if the Laplace approximation becomes more biased as rho and block size grow, it could differentially affect joint and pointwise scores and create a spurious joint advantage. The paper itself states that the Laplace approximation 'can be a crude posterior approximation' (Section 3, page 8). Please add a validation study comparing Laplace-based CV score contributions (both joint and pointwise) with MCMC- or exact-based contributions for a subset of SAR settings, especially at high rho and large test block sizes, and report the magnitude of any discrepancies.","section":"Section 3, Appendix D (Eqs. D.1–D.2, Figures D.1–D.2)"},{"comment":"The central empirical claim that joint scoring improves selection accuracy and reduces variability is supported only by visual comparison of point estimates, without confidence intervals or formal statistical tests. For example, in the covariate selection experiment (Section 3.1.1, Figure 2), correctness percentages over 1,000 replications are plotted with no error bars, making it difficult to judge whether the joint-vs-pointwise differences are beyond sampling noise, particularly at moderate rho values where the differences appear small. In the kernel experiments (Section 3.2, Figure 4), the paper reports 98% trimmed means and variances but again no uncertainty quantification. Please provide standard errors or bootstrap confidence intervals for the correctness percentages and Z-ratio estimates, or conduct a paired comparison across the common data draws, so the reader can assess the size and precision of the reported effects.","section":"Section 3, Figures 2–4"}],"minor_comments":[{"comment":"The text below Eq. (6) refers to 'the sample standard error and mean' as defining the Z ratio, but the denominator of Eq. (6) is the sample standard deviation, not the standard error; please correct the terminology.","section":"Section 3.1.1, Eq. (6) and surrounding text"},{"comment":"The caption's panel references appear to be mismatched: the text refers to 'Panel (1a)' and 'Panel (1b)' for both correctness and Z-ratio panels, and the parenthetical descriptions do not match the figure layout; please revise for clarity.","section":"Figure 4 caption"},{"comment":"The conclusion states that leave-group-out CV with 'a dozen or few dozen folds' performs well, but the simulations include test set sizes up to ntest = 144, which corresponds to only 4 folds; please qualify the claim or exclude the very large block sizes.","section":"Section 3.2 and Section 5"},{"comment":"The threshold 'standardized rho > 3/4' for a clear joint advantage appears stronger than what the simulations show, since joint evaluation also appears to help for rho = 0.5 in the network selection experiment (Figure 3) and in some kernel settings (Figure 4); please align the wording with the full range of results.","section":"Abstract, Section 3.1.2, and Section 5"},{"comment":"The definition of \\hat Z_y uses within-fold variance from a single dataset, whereas the Z ratio in Section 3 is computed across independent data draws; this distinction is not explained, so please clarify the relationship between the two quantities.","section":"Section 4, Eq. (10)"},{"comment":"The paper would benefit from a statement on code and data availability, since the simulation study is computationally intensive and reproducibility is otherwise difficult to assess; also, there is a minor typo in Section 5 ('invovling' should be 'involving').","section":"Section 5 and Appendix D"}],"recommendation":"major_revision","confidential_remarks":"The paper is well-written and the research question is timely. My main technical concern is that the headline claim rests entirely on simulation results computed with a Laplace approximation that is validated only for a single draw of kernel models, not for the SAR experiments or the held-out predictive densities; the absence of uncertainty quantification in the simulation comparisons compounds the problem. These issues are addressable within the paper's scope by adding targeted MCMC validation for the SAR settings and standard errors or confidence intervals for the reported estimates. If the authors can provide that evidence, the paper could be suitable for publication in a leading statistics journal; otherwise the claims should be substantially softened."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"The take-home: this is a solid, modest extension of Cooper et al. (2024) from time series to spatial block CV, and the main recommendation—use joint scoring when spatial dependence is strong—is plausible and supported by several independent experiments. It is not a breakthrough, but it is a genuinely useful workflow result for applied spatial statisticians.\n\nWhat is new and good: the paper systematically compares pointwise versus joint evaluation of the log score in leave-group-out CV for SAR covariate selection, network structure selection, and covariance kernel selection. The simulation design is clear, with varying test-set sizes and standardized rho values, and the Z-ratio summary is a sensible way to separate location from variability. The case studies with real data (child vaccination, lung cancer) are a nice addition; they use INLA with actual posterior draws rather than the Laplace shortcut, and they agree with the simulations. The authors are also honest about limitations: no iid comparisons, only the log score, and the difficulty of generalizing beyond standardized SAR rho.\n\nThe soft spots: the main one is the Laplace approximation. All of the simulation-based CV scores are computed via the Laplace approximation described in Appendix D, but the only validation of that approximation is for two kernel-model marginal posteriors on a single data draw. There is no check for the SAR experiments that deliver the headline rho > 0.75 result, and no check on the held-out predictive densities p(ytest | ytrain) that actually feed the joint scores. Since joint scores use the full predictive covariance while pointwise scores use only marginals, an approximation error in the covariance could affect the two scoring rules asymmetrically—and that asymmetry could in principle create the very effect the paper reports. I am not saying that is what happened; the case studies, which do not use Laplace, point the same direction. But the simulation evidence would be much more convincing with a few MCMC checks for the SAR predictive densities, especially at high rho and larger test blocks. That is a real gap, not a technicality.\n\nMinor points: the covariance kernel experiments (Section 3.2) are noisy and use trimmed means without formal statistical comparison, but the SAR experiments are clean and carry the main argument. No code is shipped, which makes the results harder to reproduce and check.\n\nThe paper deserves serious peer review. I would recommend conditional acceptance after the authors add at least a spot-check of Laplace against MCMC for the SAR setting and for the joint predictive densities, and ideally make the simulation code available. If those checks confirm the simulations, this will be a cite-worthy contribution to the spatial CV workflow literature.","headline":"Useful, honest extension of joint CV scoring to spatial models, but the simulation evidence leans on a Laplace approximation that gets only thin validation.","tokens_in":27030,"tokens_out":1477,"would_cite":true,"duration_ms":15987,"reading_group":"yes","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":["62M30","62F15","62-08"],"pacs":[],"model":"deepseek-v4-flash","headline":"Scoring left-out spatial blocks jointly, rather than pointwise, makes cross-validation model selection more accurate when spatial dependence is strong.","keywords":["spatial cross-validation","joint scoring","leave-group-out cross-validation","Bayesian model selection","simultaneous autoregressive model","Gaussian Markov random field","block cross-validation","predictive scoring rules"],"falsifier":"Re-run the pointwise-versus-joint comparison on a strongly dependent spatial process outside the Gaussian SAR family with row-standardized adjacency, such as a CAR or BYM2 model on an irregular geography: if joint scoring does not reduce the cross-validation statistic's variance or improve correct-selection rates, the claimed effect is specific to the SAR framing rather than spatial dependence itself.","tokens_in":26152,"feed_emoji":"🗺️","tokens_out":5424,"duration_ms":48034,"temperature":0.7,"pith_summary":"This paper argues that cross-validation model selection for spatial data should score each left-out block as one multivariate prediction rather than as a collection of pointwise predictions. Simulations with Gaussian spatial models show that joint scoring reduces the variability of the cross-validation selection statistic and improves the probability of selecting the correct model under strong spatial dependence, roughly when standardized spatial dependence exceeds $3/4$. The finding matters because spatial cross-validation is a common tool for choosing among models, yet the scoring rule has received less attention than the blocking design.","feed_headline":"Joint scoring beats pointwise CV for spatial models","feed_subtitle":"Block-wise joint evaluation cuts variability and improves model choice once spatial dependence is strong.","key_machinery":"The central object is the joint predictive density $p(y_{\\text{test}_k} \\mid y_{\\text{train}_k}, M)$ for each left-out block, compared with the pointwise aggregation $\\sum_i \\log p(y_{\\text{test}_{k,i}} \\mid y_{\\text{train}_k}, M)$. The paper summarizes model-selection performance through the Z ratio, the mean of the cross-validation selection statistic divided by its standard deviation across independent data draws; joint scoring raises this ratio by reducing variance.","core_discovery":"The paper claims that jointly evaluated scoring rules outperform pointwise aggregation in block-wise spatial cross-validation. For Gaussian SAR and covariance-kernel models on a regular lattice, evaluating the joint predictive density of each left-out block yields higher correct-selection rates and higher ratios of mean to standard deviation of the selection statistic than summing pointwise predictive densities, once standardized spatial dependence is above $3/4$. The improvement comes primarily from lower variability of the selection statistic rather than from a large shift in its mean.","pith_inferences":["If the mechanism is not specific to the Gaussian SAR family, joint scoring could also improve spatial model selection for CAR, BYM, and Gaussian-process models whenever the joint predictive distribution is computable.","Because the paper only studies the logarithmic score, other proper multivariate scoring rules such as the energy score or variogram score could show similar or different patterns; this is a direct testable extension.","The threshold of standardized dependence $\\rho > 3/4$ is framed in the row-standardized SAR parameterization, so transferring the threshold to other model classes would require a model-independent measure of spatial dependence.","Block cross-validation with joint scoring may serve as a computationally cheaper alternative to leave-one-out CV, since fewer model fits are needed while selection power improves."],"forward_implications":["Practitioners using block cross-validation with a dozen or so folds and joint scoring can expect more reliable model selection when spatial dependence is strong.","Pointwise evaluation can lose accuracy as test blocks grow under strong spatial dependence, while joint evaluation does not show that decline.","When spatial dependence is weak, the choice of scoring method and test-block size matters little for selection accuracy.","The benefit of joint scoring is largest in marginal cases where candidate models differ only subtly in predictive performance.","Real-data examples show joint and pointwise scoring select the same preferred model, but joint scoring yields larger Z ratios when the fitted spatial dependence is strong."],"supporting_citations":[{"why":"Supplies the time-series result and simulation methodology that this paper extends to spatial dependence.","marker":"Cooper et al. (2024)"},{"why":"Defines block cross-validation with halo designs, the CV scheme used throughout the experiments.","marker":"Roberts et al. (2017)"},{"why":"Establishes the Z-ratio and Gaussian approximation for cross-validation model comparison statistics.","marker":"Sivula, Magnusson, Matamoros, et al. (2022)"},{"why":"Provides the row-standardized SAR parameterization used to define and standardize the spatial dependence parameter.","marker":"Ver Hoef, Peterson, et al. (2018)"},{"why":"Assesses spatial cross-validation blocking variants and motivates the spatially clustered folds used in the case studies.","marker":"Mahoney et al. (2023)"},{"why":"Defines the Gaussian Markov random field class that includes the SAR and covariance models studied here.","marker":"Rue and Held (2005)"},{"why":"Formalizes proper scoring rules, justifying the logarithmic score used for joint and pointwise evaluation.","marker":"Gneiting and Raftery (2007)"},{"why":"Provides the pointwise elpd and leave-one-out CV baseline that joint evaluation is compared against.","marker":"Vehtari, Gelman, and Gabry (2017)"}],"fun_headline_variants":["Block-wise joint scoring stabilizes spatial CV","Joint scoring cuts CV variability in spatial models","Spatial CV: joint scoring beats pointwise aggregation","Joint leave-group-out CV improves spatial model choice","Joint scoring wins for strong spatial dependence"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The evidence rests on Gaussian spatial models with a known row-standardized adjacency structure, and assumes the halo-and-block design keeps training and test sets nearly independent so the Gaussian approximation to the selection statistic is valid.","fun_headline_variants_meta":{"raw":{"variants":["Block-wise joint scoring stabilizes spatial CV","Joint scoring cuts CV variability in spatial models","Spatial CV: joint scoring beats pointwise aggregation","Joint leave-group-out CV improves spatial model choice","Joint scoring wins for strong spatial dependence"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000565,"raw_usage":{"total_tokens":2600,"prompt_tokens":790,"completion_tokens":1810,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":406,"completion_tokens_details":{"reasoning_tokens":1741}},"tokens_in":406,"tokens_out":1810,"duration_ms":11102,"temperature":1.0,"reasoning_tokens":1741,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-16T11:21:58.090735+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Re-run the pointwise-versus-joint comparison on a strongly dependent spatial process outside the Gaussian SAR family with row-standardized adjacency, such as a CAR or BYM2 model on an irregular geography: if joint scoring does not reduce the cross-validation statistic's variance or improve correct-selection rates, the claimed effect is specific to the SAR framing rather than spatial dependence itself.","supporting_citations":[],"review_version":1}