Corrections to: Changing Climates of Conflict: A Social Network Experiment in 56 Schools

The problem

The article says that the school-level confidence intervals and $P$ values use “robust SEs under a normal approximation”, and it does not name the estimator. The authors' script 4. School results-with-correction.R (lines 162–163 and 208–209) uses the heteroskedasticity-robust variance with no finite-sample adjustment, known as HC0. The wild bootstrap of Tables S6 and S7 (lines 178–195 and 228–248) multiplies the unscaled residuals by standard normal draws. Given the data, its variance is exactly HC0. The row “Wild Bootstrap SE” is therefore not a second check on the row “Heteroskedasticity Robust SE”.

The regressions for disciplinary reports have 49 schools and many coefficients: block effects for 14 blocks, treatment, the composition terms and covariates. On the public data the regression of Table S6 has 19 estimated coefficients and that of Table S7 has 20, which leaves 30 and 29 residual degrees of freedom. Under a linear model with independent mean-zero errors of equal variance, fitting shrinks the residuals: their expected mean square is $30/49$ and $29/49$ of the error variance, whereas the expected HC0 variance is 0.62 of the variance of the treatment coefficient and 0.56 of the variance of the interaction coefficient. A normal reference does not correct for this, and the test rejects too often in the null experiments below.

Both results on disciplinary reports rest on this test: the main effect ($P<0.05$) and the interaction with the proportion of social referent seeds ($P=0.07$). In null experiments on the design reconstructed from the public data, the test rejects a true null hypothesis at the nominal 5% level in 12 to 16% of draws. The two printed $P$ values are weaker evidence than the same values from a calibrated test. On the public data, the test of the article gives $P=0.099$ for the main effect and $P=0.069$ for the interaction, and about $0.20$ and $0.15$ when the schools are re-randomized under the hypothesis of no effect on any school.

The description in the Methods is accurate as far as it goes, and it claims no exact size. Whether the printed $P$ values are the values of this test on the authors' own files remains unknown, because those files are not public. No replacement test is established as calibrated for the two hypotheses of the article.

Evidence a derivation

I checked the authors' script, derived the covariance identities below, calculated the estimates with the script on the public release of the data (ICPSR 37070, version 2), and examined two null experiments on the school-level design reconstructed from that release.

The bootstrap reproduces HC0. Let $X$ be the matrix of regressors, of full column rank, and write $A=(X^{\top}X)^{-1}X^{\top}$, $H=XA$ and $e=(I-H)y$ for the residuals. The bootstrap of the script regresses $X\widehat b+\operatorname{diag}(e)v$ on $X$, where the entries of $v$ are independent with mean 0 and variance 1. Its coefficient is $b^{*}=\widehat b+A\operatorname{diag}(e)v$, and hence

HC0 is too small in this design. Suppose that $y=Xb+u$, where, conditional on $X$, the errors $u$ are independent with mean zero and common variance $\sigma^2$, and let $h_{ii}$ be the $i$th diagonal entry of $H$. Then $\mathbb{E}[e_i^2\mid X]=\sigma^2(1-h_{ii})$, so

For the $j$th coefficient, the expected HC0 variance divided by the true variance is

On the public data this ratio is $0.619$ for the treatment coefficient of the regression of Table S6, where the mean of $h_{ii}$ is $19/49$, and $0.556$ for the interaction coefficient of the regression of Table S7, where it is $20/49$.

The estimates on the public data. The public-data reconstruction applies the authors' rule for selecting records, builds the network and the outcome as their scripts do, and fits their two covariate-adjusted regressions. The rule selects 49 schools, 26 control and 23 treatment, as in the article (the rule is the subject of P037-N007). The main effect is $-0.0573$ and the interaction coefficient is $-0.841$. The script prints estimates, standard errors and intervals, and no $P$ value; I computed each two-sided $P$ value here from the estimate and its standard error.

The null experiments. The first experiment draws independent standard normal errors as the outcome, 100,000 times for each of the two fixed designs, so that every coefficient is truly zero. The test of the article rejected at the nominal 5% level in 13.7% of draws for the main effect and in 16.0% for the interaction (Monte Carlo standard errors of about 0.11 percentage points). The second experiment assigns two of the four schools of each block to treatment again, 30,000 times over the 56 schools, holding the outcomes, the covariates, the draws of the seeds and the set of 49 schools fixed, as they would be if treatment had no effect on any school. The test rejected in 13.0% and 12.4% of assignments (Monte Carlo standard errors of about 0.19 percentage points). The last row of the table estimates the share of assignments in which the HC0 statistic was at least as large in absolute value as the observed one, using $(1+\text{count})/(30000+1)$.

In the first experiment the variants also miss the nominal level. HC1 with a $t$ reference rejects in 5.2% and 6.2% of draws, HC2 with a $t$ reference in 5.4% and 5.7%, and HC3 with a normal reference in 2.0% and 1.9%. The experiment has normal errors of equal variance, and it does not show how these variants behave for the randomization of the study. The re-randomization calculation is for the hypothesis of no effect on any school, with the outcomes, the covariates, the draws of the seeds and the set of 49 schools held fixed. That hypothesis is stronger than a zero average effect, and stronger than a zero interaction when treatment has some effect.

What the printed tables alone allow. An effect of $-0.0649$ with HC0 standard error $0.0251$ rounds to the entries $-0.06$ and $0.03$ of Table S6. After the rescaling of the standard error by $\sqrt{49/(49-k)}$ (HC1), its two-sided normal $P$ value is $0.043$ with $k=19$ coefficients and $0.047$ with $k=20$, the two numbers of coefficients on the public data. The rescaling alone does not remove significance for every pair of values that rounds to the printed entries. With the rescaling and a $t$ reference with $49-k$ degrees of freedom, every such pair has $P>0.05$; the smallest values are $0.051$ and $0.055$.

Influence-score terminology and current use. The raw empirical OLS influence covariance is exactly HC0. With $n$ schools, design row $x_i^{\top}$ and $\widehat\phi_i=(X^{\top}X/n)^{-1}x_i e_i=nA_{\cdot i}e_i$,

An identified current DoubleML implementation also uses an unadjusted empirical influence-score second moment: its nonclustered variance helper computes $\operatorname{mean}(\psi_i^2)/(N\operatorname{mean}(\psi'_i)^2)$, and its official variance guide uses normal inference under regularity conditions. The sources are the captured stable guide and the pinned Python source. They establish documented continued use, without measuring its prevalence or identifying the stable documentation with a particular numeric package release. Cross-fitted DoubleML scores generally differ from these in-sample OLS residuals; their shared second-moment form transfers no finite-sample calibration to this school design.

The corrected passage below retains HC0. Reporting its public-data main-effect standard error as $0.035$ rather than $0.03$ gives an extra rounding digit; it does not provide a more accurate estimator of sampling uncertainty.

Limits. The values above are values on the public data. The authors' own input files are not public, and ICPSR states that some variables of the release were changed. On the public data the coefficient and standard error of Table S6 agree with the printed $-0.06$ and $0.03$, and the standard error and $P$ value of the interaction round to the printed $0.46$ and $0.07$. The interaction coefficient is $-0.841$, which rounds to $-0.84$ where Table S7 prints $-0.85$. The rejection rates are properties of the two experiments on the public-data design. They are not the rejection probability of the study under its own unknown outcomes, and the coverage of the bands of Fig. 2 was not examined.

What survives

If the variance estimator or the reference distribution is changed, the least-squares estimates, their direction, the assignment of schools and the measurement of the networks do not change. The analysis of peer-to-peer influence computes its variance from simulated assignments and does not use this test. The correction of 30% to 25% is a separate matter (P037-N001).

The evidence for the two results on disciplinary reports is weaker than printed. On the public data the main effect is not significant at the 0.05 level under the test of the article ($P=0.099$), in agreement with the last sentence of the correction of 2018; the disagreement between that sentence and the article’s $P<0.05$ is the subject of P037-N008. The interaction keeps $P=0.07$ under a test that rejects a true null hypothesis in 12 to 16% of draws. The treatment bands of Fig. 2 use the same robust variance (their construction is the subject of P037-N026). No printed claim is shown to be false on the authors' own files. No test here is established as calibrated for a zero average effect and for a zero interaction in this design; the variants in the table weaken the evidence and are not such a test.

The problem

The Bayesian bootstrap assigns random weights to the units; the weights sum to one and have the Dirichlet$(1,\ldots,1)$ law. Dirichlet weights are obtained from independent draws of the standard exponential distribution, each divided by the sum of the draws, or, for $n$ schools, as the $n$ gaps formed by $n-1$ ordered independent $U(0,1)$ draws together with the endpoints $0$ and $1$. The notes to Tables S6 and S7 describe a different rule, a weight for each school drawn from the standard uniform distribution, and the released code computes that rule. In each of 5,000 replicates it draws one independent $U(0,1)$ weight for each school, passes the weights to the regression unchanged, and keeps the coefficient; the entry is the standard deviation of the 5,000 coefficients (4. School results-with-correction.R, lines 165–175 for Table S6 and 212–226 for Table S7). The code has no ordering, no gaps and no logarithm. Uniform weights do not have the Dirichlet law, and they do not acquire it when they are divided by their sum. The notes describe the computation correctly, and the name of the rows is wrong. The four entries of Table S6 and the eight entries of Table S7 are dispersions under uniform weights, and they are not standard errors from a Bayesian bootstrap. In the comparison using 50,000 paired refits per regression on the public data, each of the twelve uniform-weight dispersions is smaller than its counterpart with Dirichlet weights, by 27% to 32%. Read as Bayesian-bootstrap standard errors, the rows understate the dispersion under that method on the public data, and they cannot serve as a check of the robust standard errors by that method.

Evidence an exact counterexample

The evidence is an exact counterexample, a reading of the released code and a recalculation on the public data. Let $U_1,\ldots,U_{56}$ be independent and uniform on $(0,1)$, and set $P_i=U_i/\sum_j U_j$. Let $E_1,\ldots,E_{56}$ be independent standard exponential variables, and set $B_i=E_i/\sum_j E_j$, so that $(B_1,\ldots,B_{56})$ has the Dirichlet$(1,\ldots,1)$ law. Consider the event $p_1/(p_1+p_2)<1/4$. Since $P_1/(P_1+P_2)=U_1/(U_1+U_2)$, the event holds for $P$ exactly when $U_1<U_2/3$, and I have

For $B$ the event is $E_1<E_2/3$, and a direct computation gives

Equal laws assign equal probabilities to every event, and the two laws are therefore different. Only two of the weights enter the event, and the conclusion holds for any number of schools from two upwards, including the 49 schools of the peer-conflict regressions.

Normalization does not remove the difference. For a design matrix $X$ of full column rank, positive weights $w$ with $W=\operatorname{diag}(w)$, and a constant $c>0$, the weighted least-squares coefficients satisfy

A refit with uniform weights has the same coefficients whether or not the weights are divided by their sum.

For independent weights, to first order the covariance of the refitted coefficients is $\sigma^2/\mu^2$ times the unadjusted heteroskedasticity-robust covariance (HC0) of the least-squares coefficients, where $\mu$ and $\sigma^2$ denote the mean and variance of a weight. The factor is $1/3$ for $U(0,1)$ weights and $1$ for standard exponential weights. This is a statement about the linearized coefficients. It is not a correction for the sample at hand, and it does not justify multiplying the printed entries by $\sqrt{3}$.

I ran the two loops of the released code in a reconstruction on the public data (ICPSR 37070 V2), with adapted inputs and 5,000 replicates as in the code, and then refitted each regression 50,000 times with uniform weights and with Dirichlet weights built from the same draws ($U$ for the uniform weight and $-\log U$ for the exponential one). The table reports the results, and the figure shows the ratios.

The twelve ratios lie between 0.683 and 0.730, so the standard deviations under uniform weights are 27.0% to 31.7% smaller. The public data reproduce the printed tables closely but not exactly. Nine of the twelve point estimates and eleven of the twelve robust standard errors round to the printed values, and each of the twelve standard deviations under uniform weights is compatible with the printed entry once rounding and two Monte Carlo standard errors are allowed (eight of them with rounding alone). The robust covariance of the code is HC0.

What survives

This finding does not require changing the point estimates of Tables S6 and S7, the rows of heteroskedasticity-robust standard errors or the rows of standard errors from the wild bootstrap. For school-level effects, the main text computes confidence intervals and $P$ values from robust standard errors (p. 568), and the intervals of Fig. 2 use heteroskedasticity-robust standard errors (caption, p. 569); neither uses the rows in question. The student-level Wald tests use randomization-based variance estimates (p. 569); those analyses and the headline estimates do not depend on these rows. What changes is the reading of the twelve entries, which cannot be cited as standard errors from a Bayesian bootstrap. The values with Dirichlet weights in the table are Monte Carlo estimates on the public data. They are not replacements for the printed entries, which would require the authors' own data file; that file is not public, and ICPSR states that some variables of the release were changed. The comparison is specific to these regressions. It establishes no general ordering of the two dispersions and no coverage property of either weighting under the blocked design.

The problem

The article reports the effect on disciplinary reports as a 30% reduction in four places. The figure is the quotient $0.06/0.20$ of the two numbers in the Results. According to the correction of 2018, the first is the covariate-adjusted estimate of the effect and the second is the unadjusted mean of the control schools, and both are rounded. The numerator and the denominator are on different bases. The correction states that the same quotient is 29% without rounding, that the reduction is 25% when the effect and the control mean are both covariate-adjusted, and that it is 23% when both are unadjusted. From the estimated totals in the Discussion, I also obtain 25% after rounding. A citation of the original figure overstates the reduction by five percentage points, one sixth of the printed 30%.

Evidence a recalculation from the source

The evidence is a recalculation from the published numbers, a check of the pages against the correction of 2018 (PNAS 115(15):E3598), and a calculation with the authors' analysis scripts on the public release of the data (ICPSR 37070, version 2).

The published numbers. The Discussion estimates 2,695 disciplinary events without the intervention and 2,012 with it, for the same 11,938 students. The reduction is

This lies between $24.5\%$ and $25.5\%$, and it rounds to $25\%$. It equals $30\%$ only if $10\cdot 683=3\cdot 2695$, that is, $6830=8085$. If each total was rounded to a whole number of events, the reduction lies between $1364/5389\approx 25.31\%$ and $152/599\approx 25.38\%$. The quotient $0.06/0.20$ equals $30\%$ exactly, and the printed figure is not a slip of division.

The scripts and the public data. The authors' script 4. School results-with-correction.R computes the percentages of the correction. Line 208 fits the covariate-adjusted regression with the interaction of treatment and seed-group composition. Lines 253–256 average its predictions under control and under treatment over the 49 schools of the analysis, and line 286 divides their difference by the average prediction under control. I ran these lines on the public data, with an adapter that maps its variables to the inputs the scripts expect. The average predictions are $0.222537$ under control and $0.165975$ under treatment, and

The same calculation gives $28.68\%$ for the quotient of the adjusted difference and the raw control mean $0.197241$ (line 270), and $22.87\%$ for the quotient of the raw difference and the raw control mean (line 289). These round to the 25%, 29% and 23% of the correction. The raw means use all 56 schools, as discussed in P037-N009.

The regression of Table S6 is a different fit (line 162, without the interaction). On the public data its coefficient is $-0.057272$, and the difference of the average predictions of the interaction regression is $-0.056562$. Both round to the printed $-0.06$. On the public data, I reproduce the Discussion’s totals from the regression of Table S6: its average predictions over the 49 schools, multiplied by 11,938, are $2695.4$ and $2011.7$, a reduction of $25.37\%$. The three reductions $25.34\%$ (printed totals), $25.37\%$ and $25.42\%$ are different quantities, and each rounds to $25\%$.

On the public data, the authors' rule for selecting records keeps 24,191 students in 56 schools, 49 schools enter the regressions of disciplinary reports, and the treatment schools have 11,938 students. The recomputed quantities agree with the corresponding published values at their printed precision. The authors' own input files are not public, and ICPSR states that some variables of the release were changed. The values above are values on the public data. They are not the unrounded values of the authors' files.

What survives

This correction leaves the estimate of 0.06 fewer disciplinary reports per student (Table S6: $-0.06$, standard error $0.03$), the totals in the Discussion, the direction of the effect and the other analyses unchanged. The reduction is five percentage points smaller than printed. The correction requires no new causal assumption or estimator. The correction of 2018 also contains a sentence on statistical significance, which is the subject of P037-N008; the school-level tests are the subject of P037-N006.

The problem

The article says that the analyses of disciplinary reports of peer conflict (the conflict analyses) use the 49 schools for which administrative data were available, and the SI Appendix (p. 10) says that seven of the 56 schools provided no data. The authors' code does not select schools in this way. It forms each student’s peer-conflict count as the sum of 17 disciplinary categories, with a missing category counted as zero, and it keeps a school when the mean of this count over the school’s students is greater than $0.02$. This rule is a condition on the outcome. On the public data (ICPSR 37070 V2) the rule selects 49 schools, 26 control and 23 treatment, and the printed count is not shown to be wrong. The same data have disciplinary records for 52 schools, 26 in each arm. Four schools, two in each arm, have no disciplinary entry. The three further schools that the rule removes are all treatment schools. Each has a recorded count for every student, and their school means, $0.0018$, $0.0145$ and $0.0199$, are positive and below $0.02$. On the public data, “provided administrative data” does not describe the selection.

There are two consequences. First, a rule that keeps a school only when its outcome is large enough can depend on the assignment, since the treatment may change the count. A sufficient additional condition for interpreting the conflict estimates causally is that school inclusion and reporting remain unchanged across assignments, together with the paper’s model conditions; an analysis that accounts for selection could also justify that interpretation. The observed data cannot check inclusion under other assignments. Second, the estimates depend on the rule. When the three schools are admitted on the public data, the main effect moves from $-0.057$ to $-0.071$, and the interaction with the proportion of social referents moves from $-0.84$ ($P = 0.07$) to $-0.47$ ($P = 0.26$).

Evidence a recalculation from the source

The evidence is a reading of the authors' code and a recalculation on the public data. In the cleaning script (2. Cleaning.R, lines 35–38), the count is the sum of the 17 categories after missing values are set to zero; the list of categories is the subject of P037-N018. The school-level script (4. School results-with-correction.R, lines 151–153) keeps the schools whose mean count exceeds $0.02$, and it fits the regressions of Tables S6 and S7 on them. The network script (3. Network results.R, lines 167–173) applies the same rule. The script for the check within treatment schools (Fig. S5) removes the schools with mean below $0.02$, which selects the same schools here, since no school mean equals $0.02$. The rule is not applied to the survey outcomes.

I ran these lines on the public data, after converting its variable names and codes to the form that the scripts read. The authors' filter on the records keeps 24,191 students in 56 schools, as published.

The counts of the three treatment schools that have records and fail the rule sum to 13. Two of these schools have ordinary amounts of other disciplinary records: $0.51$ and $0.41$ per student, which exceed the rates of 32 and of 28 of the 49 retained schools. Nearly all of these records (245 of 255, and all 228) are generic rule violations, and the third school has one other record. These totals do not show whether the three schools reported peer conflict completely. The public data also do not show which seven schools originally provided no data.

I then fitted the covariate-adjusted regressions of Tables S6 and S7 on the public data twice, changing only the schools included.

These $P$ values use the article’s HC0 normal test. In null experiments on the 49-school public-data design, that test rejects a true null too often (P037-N006). The row for 52 schools is a check of sensitivity and not a corrected estimate.

A constructed example shows that the further condition does not follow from randomization, even when every record is complete. Suppose that a block of four schools has control means $0$, $0.03$, $0.06$ and $0.09$ per student and a constant effect of $+0.03$, and that two of the four schools are treated at random. The first school passes the rule only when it is treated. Over the six assignments, the difference in means between the retained treatment and control schools is $-0.03$, $0$, $0.03$, $-0.015$, $0.03$ and $0.075$. The average is $0.015$, half the effect.

The public data reproduce the published 24,191 students, the 49 schools and their division into 26 and 23, and the $-0.06$ and $0.03$ of Table S6. For the interaction of Table S7, I find $-0.84$ against the printed $-0.85$, with the printed standard error $0.46$. The authors' own data file is not public, and ICPSR states that some variables of the public data were changed. The values above are values for the public data, and not the values that the authors would obtain from their own files.

What survives

The number 49 and its division into 26 control and 23 treatment schools are recovered by the authors' rule on the public data. The assignment of the 56 schools is unaffected, the rule is not applied to the survey outcomes (talking about conflict, norms, wristbands), and the main effect on disciplinary reports is negative under both rules. No published estimate is shown to be wrong. What changes is the description of the sample in every analysis of disciplinary reports (Figs. 1, 2D and 3D, the peer-conflict columns of Tables S5–S7, and the checks of the SI Appendix that use the same rule). Interpreting the conflict estimates causally needs the additional justification stated above. Whether seven schools originally provided no data, and whether that failure was unrelated to treatment as the SI Appendix says, cannot be decided without the original records of the schools.

The problem

Both sentences describe the overall school-level effect of the intervention on disciplinary reports of peer conflict. The article reports this effect as a significant reduction, with $P < 0.05$ from a heteroskedasticity-robust Wald test. The correction changes the size of the reduction from 30% to 25% and ends by saying that the overall effects remain insignificant at the $0.05$ level. It names no new test, and the four passages that it lists for amendment do not include the sentence that reports significance. For one test, the two statements cannot both hold. The word “remain” is also inaccurate, since the article did not report the overall effect as insignificant. The corrected publication still carries “significantly reduced ($P < 0.05$)” in the Results beside a correction that says the opposite, and the abstract, the Significance statement, the Discussion and the legend of Fig. 2 can still be read as establishing the overall reduction at the $0.05$ level.

The authors' code identifies the analysis; its reconstruction on the public data (ICPSR 37070 V2) supports nonsignificance under a conventional two-sided Wald test. The code fits the regression of Table S6 on the 49 schools of the conflict analyses and prints the coefficient of treatment with its heteroskedasticity-robust (HC0) standard error. On the public data the estimate is $-0.0573$ reports per student with standard error $0.0347$, and the two-sided normal Wald $P$ value is $0.099$. This supports the correction’s nonsignificance statement for that reconstruction.

Evidence a derivation

The evidence is a check of the pages, a reading of the authors' code and a recalculation on the public data. Both quotations are as printed, and the correction specifies no replacement test. In the school-level script (4. School results-with-correction.R), lines 162–163 regress the school mean of the peer-conflict count on treatment, block fixed effects, the proportion of social referents among the seed-eligible students and three school covariates, with HC0 standard errors. The script prints no $P$ value. I ran these lines on the public data, after converting its variable names and codes to the form that the scripts read.

For peer conflict on the public data, I obtain

and the 95% interval is $[-0.125, 0.011]$. A separate calculation of the same regression by matrix algebra agrees to ten decimals. The code uses another fit for the corrected percentage, the interaction model averaged over the 49 schools; on the public data its contrast is $-0.0566$ with standard error $0.0321$ and $P = 0.078$. A one-sided $P$ value for the public-data coefficient of Table S6 would be $0.0495$, but neither the article nor the code specifies a one-sided test.

The rounded entries of Table S6 do not decide the question by themselves: an estimate of $-0.057$ with standard error $0.028$ or $0.034$ prints as $-0.06$ and $0.03$, and its $P$ value is $0.042$ or $0.094$.

The public data reproduce all eight entries of the table above at the printed precision, and the published 24,191 students, 56 schools and 49 schools of the conflict analyses. The authors' own data file is not public, and ICPSR states that some variables of the public data were changed. The $P$ value of $0.099$ is the value for the public data, and the $P$ value of the original run is not recovered.

What survives

The point estimate ($-0.06$ reports per student, standard error $0.03$) is reproduced on the public data, and so is the corrected 25% ($25.4\%$). Nothing here bears on the randomization, the intervention or the data. What changes is the sentence of the Results; readings of the abstract, the Significance statement, the Discussion and the legend of Fig. 2 as establishing the overall reduction at the $0.05$ level need the same qualification for the public-data analysis. Whether this test is reliable with 49 schools is the subject of P037-N006. The results on the other school-climate outcomes, on the composition of the seed group and on peer-to-peer influence rest on separate estimates and tests.

The problem

The article restricts the analyses of disciplinary reports to 49 schools (p. 567); the authors' code selects 49 schools on the public data by a rule on the school mean. The lines of the code that compute the raw means of the two arms do not apply the restriction. They average over all 56 schools. On the public data (ICPSR 37070 V2), missing disciplinary fields are coded as zero: four of the seven excluded schools have no recorded fields for peer conflict, and three have positive counts. The two excluded control schools have no recorded event, and the five excluded treatment schools have 13 events in all. Including them lowers the treatment mean more than the control mean. On the public data the means over the 56 schools are $0.197$ and $0.152$, a reduction of 22.9%, which reproduces the printed $0.20$ and 23%. Over the 49 schools the means are $0.212$ and $0.184$, a reduction of 13.6%. The correction places the 23% beside the adjusted 25% as the same comparison without adjustment. On the 49-school sample selected by the code in the public data, the raw reduction is about half of the adjusted one.

Evidence a recalculation from the source

The evidence is a reading of the authors' code and a recalculation on the public data. In the school-level script (4. School results-with-correction.R), lines 151–153 select schools whose mean exceeds $0.02$ (P037-N007), yielding 49 schools on the public data, and the fits on lines 162 and 208 use them. Lines 264, 273 and 289 compute the raw means and the raw percentage without the rule. Lines 253–256 and 286 compute the adjusted percentage from the interaction model, by averaging its predictions over the 49 schools with treatment set to zero and to one. I ran these lines on the public data, after converting its variable names and codes to the form that the scripts read. Each mean weights the schools equally.

On the public data, the model of Table S6 has fitted rates $0.2258$ and $0.1685$. Multiplying them by the 11,938 students across all 28 treatment schools, I obtain the totals 2,695 and 2,012 of the Discussion after rounding. These totals are fitted values and not raw means. The public data reproduce the published 24,191 students, the 26 control and 23 treatment schools, the $-0.06$ of Table S6 and the rounded $0.20$, 29%, 25% and 23%; the interaction is $-0.84$ against the printed $-0.85$. The authors' own data file is not public, and ICPSR states that some variables of the public data were changed. The decimals above are values for the public data; the exact historical values and excluded-school membership cannot be established without the authors' original files.

What survives

The adjusted estimate and the corrected 25% are reproduced at published precision on the public data, as are the 29% of the correction and the totals of the Discussion. No defect in the regressions follows from this sample mismatch. What changes is the interpretation of the raw control mean of the Results and the raw percentage of the correction: the code computes them over all 56 schools, rather than its 49-school conflict-analysis sample.

The problem

The released code builds the outcome as the sum of seventeen administrative counts for each student, with every missing count set to zero: the thirteen categories of the SI list, two categories that the list omits (P037-N018), the H.I.B. count as offender of item 27, and the H.I.B. count as target of item 28. The count of item 28 does not record discipline of the student. The investigators' codebook defines the variable (hibt) in these words: “Since HIB incidents were reported by offender, if 3 offenders targeted one other student as part of one event, the target would be listed in 3 writeups, and this value would reflect that number of writeups (3), not the number of incidents (1)”. A student therefore scores on the outcome each time a writeup names the student as the target of harassment by others. On the public data this count contributes 903 of the 4,075 units of the outcome (22.2%). There are 526 students with a positive count as target, and 389 of them (74%) have zero in the other sixteen counts, so that their outcome is positive through writeups as target alone. The main text describes the level of the outcome as the number of times each student was disciplined (0.20 in control schools, 0.08 at 20% social referents), and the caption of Fig. 1 describes the colours and node sizes in the same way. These descriptions are inaccurate. The outcome is an index of administrative records of peer conflict, which counts disciplinary reports of a student and writeups that name the student as target. Item 28 also misstates its own unit: the count is a number of writeups, and it equals the number of incidents represented for that student only if each such incident contributes exactly one writeup naming that student.

Evidence an exact counterexample

The evidence is a check of the codebook, an exact count, a reading of the released code and a recalculation on the public data.

For a student $t$ and an incident $e$, let $m(e)$ denote the number of offender writeups for $e$ that name $t$ as target. The variable equals $W=\sum_e m(e)$, and I write $I$ for the number of incidents with $m(e)>0$, the incidents represented in these writeups. Subtracting, I obtain

with equality only if $m(e)=1$ for every incident with $m(e)>0$. In the codebook’s example one incident has three offenders, each of whom names the same target, and I find $W=3$ and $I=1$.

The script 2. Cleaning.R, lines 36–38, forms the outcome from the seventeen variables de4 to de16, de46, de47, hibo and hibt, and the school and network scripts use this variable when the disciplinary outcome is selected. I ran these lines on the public data (ICPSR 37070 V2), with an adapter for the release’s encodings; the filter of the first script retains 24,191 of the 24,471 records, in 56 schools.

The 389 students whose outcome is positive through the count as target alone account for 650 of the 903 writeups, and all of them are in the 49 schools of the conflict analysis (26 control, 23 treatment). I then ran the school-level regressions of the script 4. School results-with-correction.R (lines 162 and 208, with the robust covariance of lines 163 and 209) on the public data, once with the outcome as coded and once without the count as target. The same 49 schools pass the script’s screen in both cases.

On the public data the main effect and its standard error round to the printed values of Table S6, and the interaction is $-0.84$ where Table S7 prints $-0.85$. The script’s reporting line (line 273) averages over all 28 control schools. Its value on the public data is 0.197, which rounds to the printed 0.20; the choice of schools in that average is the subject of P037-N009.

What survives

The survey outcomes (perceived norms, conversations with friends, wristbands), the randomization and the network definitions are not affected. The estimates for the disciplinary outcome remain estimates for a defined index of administrative records, and the printed 0.20 is reproduced to rounding on the public data as the mean of that index. Without the count as target, the two coefficients of the school-level regressions keep their sign and a similar size on the public data. This comparison does not establish that the published standard errors or statements of significance are unchanged, and the values come from the public data: the authors' own data file is not public, and ICPSR states that some variables of the release were changed. What changes is the description of the outcome and of its levels as each student’s own discipline, in the Results and in the caption of Fig. 1, and the unit of item 28.

Three points are not established. I did not recompute the network estimates of Fig. 3D and Table S5 or the drawing of Fig. 1 without the count as target. The check does not show that one incident was counted both under the code “Bully/HIB” (item 3) and under the offender count of item 27; on the public data 58 students have both counts positive, 26 only the first and 608 only the second, and these numbers do not link the counts to incidents. Removing the count as target does not turn the remaining sum into a count of distinct disciplinary actions.

The problem

The article, the SI Appendix (pp. 14, 18 and 19) and the note to Table S5 describe one population for every network result: the students with positive probability of each of the four exposures, 2,451 in number. This population is fixed before the assignment. The authors' code computes the results on another set, which I call the analysis set of the outcome. For each outcome the code keeps a student when the preprocessed pretreatment score is available and four conditions hold: (i) each of the four exposure probabilities, as estimated from $10{,}000$ simulated assignments, is greater than $0.01$; (ii) the student was not a seed in the realised assignment; (iii) the outcome is observed; and (iv) for disciplinary reports, the school’s mean count exceeds $0.02$ (P037-N007). The regression is fitted on the analysis set, and the fitted values are averaged over the same set.

Condition (i) is stricter than positive probability. Condition (ii) depends on the realised assignment, and conditions (iii) and (iv) depend on the outcome. The SI Appendix sets the directly treated exposure aside (p. 17), which accounts for leaving the outcomes of the realised seeds out of the fit. It does not say that these students are also left out of the set over which the means are averaged. The reported means are averages over a set that depends on the assignment and on the outcome, so the number of students must be reported for each outcome.

Evidence an exact counterexample

The evidence is a reading of the authors' code and a recalculation on the public data (ICPSR 37070 V2). In the network script (3. Network results.R), line 190 sets the value $0.01$, line 195 defines the analysis set by conditions (i)–(iv) and pretreatment-score availability, line 216 fits the regression on it, and lines 223–226 average the fitted values over it. For norms the code imputes missing pretreatment scores; for the other outcomes it uses a constant. The availability condition removes no student on the public data.

An exact example shows that condition (i) differs from positive probability. Suppose that a student who is not eligible has six eligible peers, one of them a social referent, in a block of 12 eligible students from which six seeds are drawn. There are $\binom{12}{6} = 924$ equally likely draws, and the student has no seed peer in exactly one of them. With a treatment probability of $1/2$ for the school, the probabilities of the control school, of no seed peer, of a non-referent seed peer only and of a referent seed peer are $1/2$, $1/1848$, $461/1848$ and $1/4$. All are positive, and the second is below $0.01$.

I recalculated these counts on the public data using the scripts, after converting its variable names and codes to the form that the scripts read. I used the random seed of the code and the sampling method of the R versions of the time (before 3.6).

With the exact probabilities in condition (i), 2,921 of the 3,120 students pass. Using exact probabilities, the cutoff removes 199 students of the positive-probability population on the public data; the further membership changes from simulation noise are the subject of P037-N021. Condition (ii) removes 61 students. Seed-eligible students who were not drawn as seeds stay in the analysis set: 207 of the 2,847 students of the wristband analysis are seed-eligible. The code holds the analysis set fixed in the $1{,}000$ simulated assignments from which the standard errors are computed. Since no student of the analysis set has the directly treated exposure as observed exposure, the imputation of potential outcomes on SI pp. 19–20 never needs a mean for that exposure.

The public data reproduce the published 24,191 students in 56 schools. None of the four analysis sets has 2,451 students. The authors' own data file and their saved draws are not public, and ICPSR states that some variables of the public data were changed. The four numbers of the original run are not recovered, and the counts in the table are counts for the public data.

What survives

Positive probability as a requirement, the four conditions of exposure and the reported network estimates are untouched. No estimate, $P$ value or conclusion of the article is shown to change. The school-level regressions (Tables S6 and S7) do not use the analysis sets. What changes is the population to which Fig. 3 and Table S5 refer, and its size for each outcome: describing the rule correctly replaces the fixed population of the article by a set that depends on the assignment and on the outcome.

The problem

Each panel of Fig. 2 comes from a school-level regression of the outcome $Y_i$ on the treatment indicator $T_i$, its product with the proportion $C_i$ of social referents among the seeds, block indicators and school covariates $z_i$ (SI equation (2), with the adjustments of SI p. 14):

After omitting redundant design columns, let $\bar z$ be the average of $z_i$ over the schools in the fit, and write $\ell_0 = (1, 0, 0, \bar z^{\top})$, $\ell_1(c) = (1, 1, c, \bar z^{\top})$ and $d(c) = \ell_1(c) - \ell_0 = (0, 1, c, 0^{\top})$. The red line is the predicted control mean $\widehat m_0 = \ell_0 \widehat b$, where $\widehat b$ is the vector of estimated coefficients, and the black line is the predicted treatment mean $\widehat m_1(c) = \ell_1(c)\widehat b$ at proportion $c$. The caption and the keys say that the dashed curves are 95% confidence intervals for these two means, computed from the robust covariance matrix $\widehat V$ of the fit. Such intervals are

The plotting code draws something else, in four respects.

1. The treatment band. For each $c$ the code refits the regression with $C_i$ replaced by $C_i - c$ and takes the robust standard error of the coefficient on $T_i$. When the refit uses the same schools and covariates as the line, that coefficient is $\widehat\tau_0 + c\widehat\beta$, the difference between the two lines, and its variance is $d(c)\widehat V d(c)^{\top}$. The band uses a standard error of the contrast, drawn around the treatment line. It treats the control mean as known.
2. The control band. The code uses $1.96\sqrt{s_0^2/56}$, where $s_0^2$ is the sample variance of the residuals of the control schools. This is not $1.96\sqrt{\ell_0 \widehat V \ell_0^{\top}}$, and the divisor is 56 although the fit has 28 control schools, or 26 for peer conflict on the public data.
3. Panel D. The line is fitted on the 49 schools of the conflict analyses (P037-N007), but the refit that gives its treatment band uses all 56 schools.
4. Fig. S6. The refit behind the treatment bands controls for the share of social referents by indegree in the seed-eligible group, where the fit of the lines controls for the share by betweenness. The bands of Fig. S6 share the two variance errors, panel (d) also uses the mismatched sample, and its lines are placed on the control mean of the indegree regression.

The bands are not wrong in one direction. On the public data, the control band is narrower than the captioned interval in every panel. At the three proportions checked (5%, 10% and 20%), the treatment band is wider than the captioned interval in ten of the twelve cases and narrower in two. The Results cite Fig. 2 for the school-level effects, and a judgement of uncertainty from the overlap or the width of these bands does not have the meaning that the caption gives it. The coefficients, their standard errors and the tests are computed elsewhere in the code and do not use these formulas.

Evidence an exact counterexample

The evidence is a derivation with an exact counterexample, a recalculation on the public data, and measurements of the printed figures from rendered pages.

Derivation. For a fixed fitting sample and covariates, since

the refit with $C_i - c$ is the same regression in other coordinates. It has the same fitted values and residuals, its coefficient on $T_i$ is $\widehat\tau_0 + c\widehat\beta$, and the robust variance of that coefficient is $d(c)\widehat V d(c)^{\top}$. The variance of the treatment mean is

which differs from $d \widehat V d^{\top}$ by the variance of the control mean and twice a covariance. Since the covariance may be negative, the drawn band may be wider or narrower than the captioned one.

Counterexample. The synthetic example has 56 schools and no covariates: 14 control schools with $Y = 0$ and 14 with $Y = 1/5$; 7 treatment schools with $(C, Y) = (0, 1/10)$, 7 with $(0, 3/10)$, 7 with $(1/5, 2/5)$ and 7 with $(1/5, 3/5)$. The least squares fit of $Y = \alpha + \tau_0 T + \beta\, T C$ is $(\widehat\alpha, \widehat\tau_0, \widehat\beta) = (1/10, 1/10, 3/2)$, every residual is $\pm 1/10$, and the robust covariance matrix is

At $c = 1/10$ the two predicted means are $1/10$ and $7/20$. Each has variance $1/2800$, and the contrast has variance $1/1400$, with covariance $-1/2800$ between the control mean and the contrast. The code’s treatment band uses $1/1400$, twice the variance of the treatment mean. Its control band uses $s_0^2/56 = (28/100)/(27 \cdot 56) = 1/5400$, which is $14/27$ of the variance of the control mean. All of these values are exact rational numbers.

Recalculation. I recalculated the authors' school regressions on the public data (ICPSR 37070 V2; 24,191 students in 56 schools after the authors' filter), using an adapter for the public release and current R packages, and compared the standard error that the plotting code uses with the robust standard error of the predicted mean from the fit of each line. The table reports both for the four panels of Fig. 2, and the figure below draws the corresponding intervals. In panel D at $c = 0.05$, the drawn standard error is $0.0374$, where the contrast from the 49-school fit of the line has standard error $0.0482$. For Fig. S6(d) at $c = 0.20$, the drawn standard error is $0.0476$ and that of the predicted mean from the betweenness fit is $0.0656$.

Measurement of the printed figure. The authors' own analysis files are not public, and ICPSR states that some variables in the release were changed. To connect the code to the printed figure, I compared the bands that the code draws on the public data with the bands as printed in the two peer-conflict panels. In the published panel D, the half-width of the control band measures $0.0289$, and the half-widths of the treatment band at $c = 0.05$, $0.10$ and $0.20$ measure $0.0744$, $0.0572$ and $0.0917$. From the code on the public data, I obtain $0.0287$ and $0.0733$, $0.0579$ and $0.0910$. In the published Fig. S6(d), the measured half-widths are $0.0302$ and $0.0758$, $0.0559$ and $0.0938$, against $0.0302$ and $0.0755$, $0.0576$ and $0.0933$ from the code. Each pair agrees within $0.002$, the tolerance of the measurement. These measurements support the code’s construction in the two peer-conflict panels. The captioned interval for the control mean of panel D would have half-width $0.0445$ on the public data.

What survives

The lines of Fig. 2 do not change. The coefficients for treatment and for its interaction with the proportion of social referents, their robust standard errors, the Wald tests and the bootstrap calculations are computed before the plotting, and none of them uses the band formulas. The peer-exposure analyses of Fig. 3 use other calculations. In Fig. S6 the slopes and the treatment coefficients stay as they are, but the plotted levels move slightly when the lines are placed on the control mean of the betweenness regression: for peer conflict on the public data that mean is $0.2171$, where the plotting code uses $0.2225$ from the indegree regression. No reading of a panel is shown to reverse. This check does not establish the coverage of any of these intervals; the behaviour of the robust standard errors with 49 schools is the subject of P037-N006. The cells of Table S7 that differ from the fit on the public data are recorded under P037-N011.