Discrete versus continuous indicators

One of the first matters confronting researchers using CR indicators is the scaling of the indicators included in the index. In perhaps the first study to use a CR index, Sameroff et al. (Reference Sameroff, Seifer, Barocas, Zax and Greenspan1987) employed discrete indicators, each scored in dichotomous fashion. Later researchers have attempted to use indicators scaled in more continuous or quantitative fashion. Each method has strengths, and each has weaknesses, all worthy of critical appraisal.

Testing hypotheses versus interrogating a model

The standard application of hypothesis testing in psychology, using a null hypothesis statistical test (NHST) approach, is to state a null hypothesis, H₀, and a mutually exclusive and exhaustive alternative hypothesis, H_A. For example, the null hypothesis might propose that all population regression coefficients in a regression equation are simultaneously zero, and the alternative hypothesis would be that at least one population regression coefficient is not zero. Data are collected, and a regression model is estimated, allowing a test of the null hypothesis. If the null hypothesis is rejected, the alternative hypothesis can be accepted, resulting in a confirmation of the hypothesis motivating the investigation. Hence, this approach has a bias in favor of confirmation, and only the theoretically uninteresting null hypothesis is tested.

A contrasting approach has been called severe testing by Mayo (Reference Mayo2018), an approach positing that disconfirmation is the path forward in science. This approach is consistent with the ideas of Popper (Reference Popper1935/1959), Meehl (Reference Meehl1990), and others that we should test our predictions, not null hypotheses, and test them as strongly or severely as possible. Under this approach, an investigator should develop a theoretical model for phenomena in a domain and then, after collecting data, test whether that model fits the data. When testing model predictions, this testing should involve severe testing of predictions or interrogating the fit of the model. If the model fits the data adequately, the model has survived a severe test and gained inductive support. Or, if the model is rejected as having poor fit to the data (i.e., model predictions are disconfirmed), something valuable has been learned – the theory was unable to account for the data, so requires modification to be in better accord with observations.

As applied to the construction of a CR index, the standard NHST approach might test whether a CR index formed as the sum of several risk indicators was a better predictor of a negative developmental outcome than was the single best indicator of risk. The contrasting, severe testing approach would test as directly and severely as possible whether equal weighting of risk indicators was acceptable or should be rejected in favor of a more complex and adequate form of weighting. If equal weighting cannot be rejected, the “equal weights” hypothesis has been interrogated and passed a severe test, supporting this simple weighting scheme.

To weight differentially or not: Interrogating a weighting sccheme

Equal versus differential weighting of risk indicators when forming a sum is a key issue. If multiple risk indicators are on the same metric, summing the indicators may serve a legitimate scientific purpose regarding equal versus differential weighting of indicators, regardless of their degree of intercorrelation. Consider an outcome variable Y and two risk indicators, X ₁ and X ₂, which are assumed to be on the same metric (e.g., standardized to have equal means and equal variances, or both 0 – 1 dichotomies). A regression equation could be written as:

(1) $${Y_i} = {B_0} + {B_1}{X_{i1}} + {B_2}{X_{i2}} + {E_i}$$

where Y _i is the score of person i on the outcome variable, X _i1 and X _i2 are scores of person i on the two risk indicators, respectively, B ₀ is the intercept, B ₁ and B ₂ are the raw score regression coefficients for the two risk indicators, respectively, and E _i represents error in predicting Y for person i. This equation would have a squared multiple correlation, or R ², that indicates the proportion of variance in Y accounted for by the weighted predictors. [Note: to ease presentation, the subscript i for person will be deleted from subsequent equations, with no loss in generality.].

Regardless of the degree of correlation between the two risk indicators, summing the two indicators would lead to the following equation:

(2) $${Y_{}} = {B_0} + {B_C}({X_1} + {X_2}) + {E^*}$$

where ${B_C}$ represents the raw score regression weight constrained to equality across the two risk indicators, E* represents the prediction error in this equation, and other symbols were defined above. Based on parameter nesting, Equation 2 is nested within Equation 1, as Equation 2 places an equality constraint on the two regression coefficients in Equation 1. Given the nesting, one could test the difference in explained variance for the two equations with the typical F-ratio for nested regression models (cf. Cohen, Cohen, West, & Aiken, Reference Cohen, Cohen, West and Aiken2003, p. 89), as:

(3) $${F_{\left( {{p_1} - {p_2}, N - {p_1}} \right)}} = {{\left( {R_{Eq1}^2 - R_{Eq2}^2} \right)/({p_1} - {p_2})} \over {(1 - R_{Eq1}^2)/(N - {p_1})}}{\rm{ }}$$

where $R_{Eq1}^2{\rm{ \hskip1pt and\hskip1pt }}R_{Eq2}^2$ are squared multiple correlations under Equation 1 and Equation 2, respectively, ${p_1}\hskip1pt{\rm{ and }\hskip2pt}{p_2}$ are the number of regression slopes estimated in the two equations, respectively, and N is sample size. The resulting F ratio has (p ₁ – p ₂) and (N – p ₁) degrees of freedom and tests the hypothesis that the two regression slopes in Equation 1 are equal. If the F ratio were larger than the critical value at a pre-specified level (e.g., α = .05), the hypothesis of equality of regression coefficients could be rejected, and differential weighting of risk indicators is justified. On the other hand, if the F ratio did not exceed the critical value, the hypothesis of equality of regression coefficients cannot be rejected, and equal weighting is appropriate.

Note that Equations 1 and 2 easily generalize to situations with more than two risk indicators, but only if all indicators are on the same metric. If one had 10 risk indicators, an initial equation would have 10 predictors ${X_1}{\rm{ through }}{X_{10}}$ each with its own regression weight, as:

(4) $$Y = {B_0} + {B_1}{X_1} + {B_2}{X_2} + \ldots + {B_{10}}{X_{10}} + E$$

where symbols were defined above, and Equation 2 would then become:

(5) $$Y = {B_0} + {B_C}({X_1} + {X_2} + \ldots + {X_{10}}) + {E^*}$$

where ${B_C}$ represents the raw score regression weight constrained to equality across all 10 risk indices, and other symbols were defined above. The proper adaptation of Equation 3 would then provide an omnibus test of equality of the regression weights across all 10 risk indicators, an F-ratio that, in this case, would have 9 and (N – 11) degrees of freedom. If the resulting F ratio were nonsignificant, the hypothesis that the regression weights were simultaneously equal could not be rejected, justifying equal weighting of all 10 indicators when forming a CR index. Of course, if the F ratio were significant, the hypothesis of equality of regression weights would be rejectable, and a more informed and complex form of weighting might be considered. If equality of all regression weights were rejected, theory or the pattern of weights when all were freely estimated might offer options for more complex, but still restricted weighting schemes.

In evaluating comparisons such as those outlined in Equations 1 through 5, sample size and power to reject the hypothesis of equality of regression weights must be considered. As sample size increases, power to reject the “equal weights” hypothesis will increase. With extremely large sample size, an “equal weights” hypothesis might be rejectable via statistical test, even if the differences in the regression weights across predictors are trivial in magnitude. Hence, some notion of practical significance of the difference must also be weighed, whether of the magnitude of the differences in the raw score regression weights or the difference in the R ² values for the equations. Conversely, small sample size and resulting low power may lead to difficulty in rejecting the hypothesis of equality of regression weights even if these weights, in truth, differ in the population, essentially committing a Type II error. But, with small sample size, it may be more appropriate to proceed in the presence of a possible Type II error rather than promote differential weighting that cannot be justified statistically.

Equality versus differential weighting of predictors in regression analysis has surfaced as an issue with regularity over time (e.g., Wilks, Reference Wilks1938; Wainer, Reference Wainer1976; Dawes, Reference Dawes1979). Many have argued that equal weights are expected to lead to little loss in predictive accuracy relative to differential weights in many situations. Regrettably, this literature is beset with conflicting claims, much too voluminous to review here. In a different, though related vein, experts recently have discussed the fungibility, or exchangeability, of coefficients in regression (e.g., Waller, Reference Waller2008) and structural equation models (e.g., Lee, MacCallum, & Browne, Reference Lee, MacCallum and Browne2018). These researchers have cautioned against having too much faith in precise, optimal regression weights. For example, in a regression equation with three or more predictors, if a small drop in explained variance is allowed, an infinite number of different sets of regression weights all produce the same R ², and many of these sets of regression weights may have little similarity to the optimal least squares estimates in the equation that maximizes R ². Thankfully, if equal weighting of predictors is justified, the problem of fungible regression weights largely vanishes, as the number of regression weight estimates is reduced and the flexibility of the equation is curtailed.

In prior work on CR indices, equal weighting of risk indicators has virtually always been utilized. Equal weights may not be optimal in a least squares sense, but may be preferred on several grounds, including parsimony, efficiency, and openness to replication in subsequent studies. However, the issue of equal versus differential weighting of risk indicators should be interrogated or subjected to severe test, as this may impact the relations between a CR index and outcomes it should predict, so has a bearing on the construct validity of the CR index.

Three guiding principles

When formulating and then interrogating a cumulative index derived from a number of risk indicators, the weighting and summing of indicators should be guided by justifiable analytic principles. Three principles are here proposed.

Background

The study by Sameroff et al. (Reference Sameroff, Seifer, Barocas, Zax and Greenspan1987) was one of the first, if not the first, to use a CR index of environment risk for low intelligence. The sample consisted of 215 mothers and children, and child WPPSI Verbal IQ at age 4 years was the outcome variable. The 10 dichotomously scored risk factors are shown in Table 1. Certain variables (e.g., ethnic status, family support) were scored in direct and unambiguous fashion, and others (e.g., occupation, education) were dichotomized at commonly used points on their respective continua. For the remaining six variables, Sameroff et al. identified cut-scores that would leave about 25% of the sample with the highest risk receiving a risk score of 1, with the remaining 75% being assigned a risk score of 0. Additional details about each of the predictors were provided by Sameroff et al.

Table 1. Ten Risk Factors Used by Sameroff et al. (Reference Sameroff, Seifer, Barocas, Zax and Greenspan1987)

Note: Tabled material is adapted from Table 2 of Sameroff et al. with minor modification.

In Table 2, correlations among the 10 risk factors and the WPPSI Verbal IQ outcome variable are shown, along with estimated means and SDs of the variables (see Supplementary Material for how information from Sameroff et al., Reference Sameroff, Seifer, Barocas, Zax and Greenspan1987, was used to develop Table 2). Inspection of Table 2 reveals that all 10 risk factors have negative correlations with WPPSI Verbal IQ, consistent with mean differences between low-risk and high-risk groups reported by Sameroff et al. The correlations in Table 2 also exhibit a strong positive manifold of correlations, with 44 of the 45 correlations among risk factors of positive valence. The single exception was the small negative correlation, r = −.03, between maternal mental health and ethnic status.

Table 2. Correlations and Descriptive Statistics for Child Verbal IQ and 10 Risk Variables from Sameroff et al. (Reference Sameroff, Seifer, Barocas, Zax and Greenspan1987)

Note: N = 215. The correlations in the table above have been re-arranged from those reported by Sameroff et al. (Reference Sameroff, Seifer, Barocas, Zax and Greenspan1987). As explained in the Appendix, three risk factor variables were reverse scored, and the correlation between WPPSI Verbal IQ and Occupation of r = −.59 as reported in text of the Sameroff et al. article was used, replacing the r = −.58 in their Table 3.

The risk factors shown in Table 1 are a varied amalgam of risk factors. The first four risk factors in Table 1 will hereinafter be called the Classic 4, ^{Footnote 1} given substantial research published over the past half century or longer on relations of these variables to child intelligence. From early work by Terman (Reference Terman1916) and Brigham (Reference Brigham1923), down through work by Broman (Reference Broman1987; Broman et al., Reference Broman, Nichols and Kennedy1975), Jensen (Reference Jensen1998) and Herrnstein and Murray (Reference Herrnstein and Murray1994), and recently Johnson, Brett, and Deary (Reference Johnson, Brett and Deary2010), researchers have studied relations of SES (occupation), education, and ethnicity with intelligence. The fourth indicator in the Classic 4 set – maternal interaction – is also based on substantial work. Yarrow and associates (e.g., Messer, Rachford, McCarthy, & Yarrow, Reference Messer, Rachford, McCarthy and Yarrow1987; Yarrow, MacTurk, Vietze, McCarthy, Klein, & McQuiston, Reference Yarrow, MacTurk, Vietze, McCarthy, Klein and McQuiston1984; Yarrow, Rubenstein, & Pedersen, Reference Yarrow, Rubenstein and Pedersen1975) reported strong relations between parental stimulation and interaction during a child’s infancy and young child problem-solving. Other work by Bradley and Caldwell and colleagues (Bradley & Caldwell, Reference Bradley and Caldwell1980; Bradley, Caldwell, & Elardo, Reference Bradley, Caldwell and Elardo1977; Elardo, Bradley, & Caldwell, Reference Elardo, Bradley and Caldwell1975) found strong relations (e.g., correlations ranging from .36 to .66) between observed mother-child interaction during the child’s infancy and child IQ at 3 years.

In contrast to the Classic 4, the remaining six risk variables, here referred to as the Modern 6, have received much less research attention with regard to predicting children’s intelligence. Certainly, father absence, high maternal anxiety, high numbers of difficult life events, and the remaining indicators likely confer risk for poorer development in general. But, these six indicators have received far less consistent attention as predictors of low child IQ as have the first four risk indicators.

Regression analyses

Because the 215 observations had complete data on all variables, the summary data in Table 2 can be used to perform regression analyses. For all multiple regression analyses reported here, I used ordinary least squares (OLS) estimation in the PROC REG program in SAS (all analysis scripts are contained in Supplementary Material). To evaluate model fit, I used the R ² for a model and differences in R ² for competing models, the F-test for differences in fit of nested models, the adjusted R ² (adj-R ², which has a penalty for model complexity), and the Schwarz Bayesian Information Criterion (BIC). ^{Footnote 2} BIC values are not on a standardized metric, but lower values indicate better model fit.

Replicating results reported by Sameroff et al

Sameroff et al. (Reference Sameroff, Seifer, Barocas, Zax and Greenspan1987) reported that the single best predictor of Child Verbal IQ was occupation, R ² = .35. They noted that including all 10 risk factors led to a much higher level of explained variance, R ² = .51. Unfortunately, Sameroff et al. did not present any additional details, such as parameter estimates, and their SEs, for this 10-predictor model, so it is not possible to evaluate their reported results further.

I first wanted to replicate results reported by Sameroff et al. (Reference Sameroff, Seifer, Barocas, Zax and Greenspan1987) to verify that data in Table 2 were sufficient to reproduce their results. Model 1, shown in Table 3, used occupation as the sole predictor of child Verbal IQ and led to an R ² = .348. Adding the nine remaining risk factors led to Model 2 in Table 3, which had R ² = .519. Both of these models replicated closely R ² values reported by Sameroff et al. Note that, in Model 2, only the Classic 4 predictors had regression weights significant at p < .05, and none of the Modern 6 met this criterion.

Table 3. Alternative Regression Models for the Sameroff et al. (Reference Sameroff, Seifer, Barocas, Zax and Greenspan1987) Data

Note: N = 215. Tabled values are raw score regression coefficients, their SEs in parentheses, and associated t-ratios. The R ² is the squared multiple correlation; the adjusted R ² is the shrunken estimate of squared multiple correlation that adjusts for model complexity. ^a The t ratio has df equal to ν ₂ for the F ratio for the equation, shown below. With 204 or more df, the critical t value at the .05 level is 1.98. ^b For the F ratio for each model, ν ₁ = numerator degrees of freedom, and ν ₂ = error (or denominator) degrees of freedom. ^c BIC is the Schwarz Bayesian Information Criterion.

Structural equation modeling

I conducted comparable analyses to those reported in Table 3 using structural modeling software (SEM) programs Mplus (Muthén & Muthén, Reference Muthén and Muthén1998-2019) and lavaan (Rosseel, Reference Rosseel2012) in R. Essentially identical results were obtained; given limitations of space, the full set of analysis scripts and descriptions of SEM results were placed in Supplementary Material.

Background

The second example is from Masarik et al. (Reference Masarik, Conger, Donnellan, Stallings, Martin, Schofield, Neppl, Scaramella, Smolen and Widaman2014), which used a CR index of genetic risk across five SNPs. Based on prior research investigating SNPs for environmental susceptibility, Masarik et al. designated the alleles as conferring susceptibility: (a) the short (s) allele of the 5-HTTLPR SNP in 5HTT; (b) the A1 allele of the Taq1A polymorphism in ANKK1 (often called DRD2 in prior work); (c) the 7R allele of exon-3 VNTR in DRD4; (d) the 10R allele of the 5’ VNTR in DAT; and (e) the Met allele of the Val158Met polymorphism in COMT. [For ease of presentation, the italicized acronyms for each allele are used in the remainder of this paper.].

Each of the five SNPs was originally assigned scores of 0, 1, or 2 based on how many of the plasticity alleles were observed. To enable more precise tests of possible forms of gene action, I created two redefined scores for each SNP, shown in Table 4. The dominant code reflects the assumption that only a single plasticity allele is needed to drive a genetic effect, and the presence of a second plasticity allele provides no additional risk or plasticity. The recessive code represents the hypothesis that no effect of the SNP would be seen unless two plasticity alleles are present. If the dominant and recessive codes are entered as separate predictors in an equation, several informative patterns might arise: (a) significant effect of the dominant code and near zero effect of the recessive code would support a conclusion that the SNP effect was dominant; (b) significant effect of the recessive code and near zero effect of the dominant code would support recessive gene action; and (c) equal regression weights for the dominant and recessive codes would support a conclusion of linear gene action.

Table 4. Dominant and Recessive Codes for an Individual SNP

The Masarik et al. (Reference Masarik, Conger, Donnellan, Stallings, Martin, Schofield, Neppl, Scaramella, Smolen and Widaman2014) investigation was of G×E interaction, and the environmental variable analyzed here was parental hostility toward a target adolescent when the adolescent was approximately 15 years of age. Parents and the adolescent engaged in videotaped interactions as they discussed several issues, including recent disagreements. Observers coded the videotaped interactions on 1-to-9 rating scales for three items – hostility, angry coercion, and antisocial behavior – of each parent toward the adolescent. Ratings were made for each parent and during each of two interaction tasks. The ratings were averaged across parents and across tasks to form a scale score for parental hostility, which could range from 1 to 9.

The outcome variable was obtained 16 years later when the target adolescents were now young adults averaging about 31 years of age. At this measurement point, each target individual was videotaped interacting with a romantic partner using similar methods as used for the prior parent-adolescent interactions. Observers coded these interactions using the same rating scales, which were again averaged into a scale for target hostility toward romantic partner.

Target participant sex (coded 0 = male, 1 = female) was used as a covariate, following Masarik et al. (Reference Masarik, Conger, Donnellan, Stallings, Martin, Schofield, Neppl, Scaramella, Smolen and Widaman2014). To allow comparisons of results across OLS regression and ML analyses using SEM software, analyses are restricted to the 281 participants who had complete data on the variables to be analyzed. For details on measures and procedures, see Masarik et al. (Reference Masarik, Conger, Donnellan, Stallings, Martin, Schofield, Neppl, Scaramella, Smolen and Widaman2014).

Preliminary considerations

At least four concerns arise with analyzing genetic SNP effects in this data set. The first issue involves potential correlations of SNPs with other variables in the data set and among SNPs. Significant correlations of gene SNPs with variables such as parental hostility lead to difficulties in interpretation, such as passive gene-environment correlation. In Table 5, the correlations among five variables are shown: target (aged 31 years) hostility toward romantic partner, parental hostility toward target (aged 15 years), a gene SNP total score (sum across the 5 SNPs), a gene Dominant score (sum of the 5 dominant SNP scores), and a gene Recessive score (sum of the 5 recessive SNP scores). The r = .280 correlation of the first two variables is a strong inter-generational (i.e., parent behavior to offspring behavior) relation, particularly given the 16 years elapsed between the two measures. The correlations of target hostility with the three SNP scores were very small and nonsignificant, as were correlations of parental hostility with the SNP scores. The lack of any substantial correlation between gene SNP scores and environment scores simplifies theoretical and empirical considerations.

Table 5. Correlations of Target and Parent Hostility with Target Gene Indexes, based on Masarik et al. (Reference Masarik, Conger, Donnellan, Stallings, Martin, Schofield, Neppl, Scaramella, Smolen and Widaman2014) Data

Note: N = 281. Target Hostility and T Host = Target Hostility toward Romantic Partner at Age 31 years, Parent Hostility and P Host = Parent Hostility toward Target Aged 15 years, T Gene Total = Target Gene Total Index, T Gene Dom = Target Gene Dominant Index, T Gene Rec = Target Gene Recessive Index. Correlations greater than ∼ |.12| are significant at p < .05, correlations greater than ∼ |.154| are significant at p < .01.

A second issue involves correlations among SNP scores, provided in Table 6. The first five variables in the table are the dominant scores for the five SNPs, and the remaining five variables are corresponding recessive scores. The five underscored coefficients in Table 6 are correlations between the dominant and recessive scores for each SNP. These correlations are structurally positive (i.e., necessarily positive) because one cannot have a score of 1 on the recessive score unless one has a 1 on the dominant score. These scores must correlate positively, and all were statistically significant (ps < .01), although not large. The other 40 correlations shown in Table 6 averaged almost precisely zero (−.003). Only one of these 40 correlations met the α = .05 criterion for significance, so the SNP scores were essentially independent.

Table 6. Correlations among Five Gene Dominant and Recessive Codes for Target Adolescent, using Masarik et al. (Reference Masarik, Conger, Donnellan, Stallings, Martin, Schofield, Neppl, Scaramella, Smolen and Widaman2014) Data

Note: N = 281. Underscored correlations must be non-zero, given the construction of the dominant and recessive codes for individual SNPs. Correlations greater than ∼ |.12| are significant at p < .05, correlations greater than ∼ |.154| are significant at p < .01. ^a In SNP names, d = dominant code, r = recessive code.

The third issue is the form analyses might take. With a large number of SNP scores available, one could analyze (a) only a gene total score; (b) 5 linear (0,1,2) SNP scores; (c) only the 5 dominant SNP scores and/or their sum; (d) only the 5 recessive SNP scores and/or their sum; (e) the dominant and recessive total scores, or (f) 10 separate dominant and recessive SNP scores. For simplicity, only analyses using the gene total score and the 10 separate dominant and recessive SNP scores are reported here, emphasizing strengths and weaknesses of each approach. Parent hostility, the key environmental predictor, was mean-centered, to simplify interpretation.

A fourth issue is whether it is reasonable to form an equally weighted sum of a set of scores, here SNP scores, that are essentially uncorrelated. Based on homogeneity reliability theory, this makes no sense at all, as no increase in reliability would accompany the summing of uncorrelated scores. But, the equally weighted sum of a set of uncorrelated indicators might possess considerable reliability in terms of stability over time. Moreover, the sum can be used to test whether the risk indicators have significantly different weights in predicting an outcome. What is needed is a test of whether this equality constraint is reasonable or whether it costs too much in predictive power, and such tests are demonstrated in the following section.

Regression analyses

Four regression models were fit using the genetic total score as the Gene variable, with parent hostility as Environment variable, to predict target hostility with romantic partner. The first model, Model 1 in Table 7, employed Female as control variable and used the E and G main effects as predictors. In Model 1, only Female and the Environment main effect were significant, and the overall model had a fairly large R ² = .1589. Adding the G×E interaction resulted in Model 2, which exhibited a moderately higher R ² = .1734. The increase in explained variance, ΔR ² = .0145, was significant, F(1, 276) = 4.87, p = .03, and the regression coefficient for the G×E interaction, B = 0.10 (SE = 0.05), p = .03, was significant and in the predicted direction.

Table 7. Alternative Regression Models for Target Hostility using Total Gene Index for Masarik et al. (Reference Masarik, Conger, Donnellan, Stallings, Martin, Schofield, Neppl, Scaramella, Smolen and Widaman2014) Data

Note: N = 281. E (Par Hostility) stands for “Environment main effect, Parental Hostility”; G CR Index stands for “Gene main effect, a CR Index formed as the sum of scores on 5 target genes”; and GxE interaction is the product of the environmental and gene main effects. Tabled values are raw score regression coefficients, their SEs in parentheses, and associated t-ratios. The R ² is the squared multiple correlation; the adjusted R ² is the shrunken estimate of the squared multiple correlation. ^a The t ratio has df equal to ν ₂ for the F ratio for the equation, shown below. ^b For the F ratio for each model, ν ₁ = numerator degrees of freedom, and ν ₂ = error (or denominator) degrees of freedom.

Belsky and colleagues (Belsky, Pluess, & Widaman, Reference Belsky, Pluess and Widaman2013; Belsky & Widaman, Reference Belsky and Widaman2018; Widaman, Helm, Castro-Schilo, Pluess, Stallings, & Belsky, Reference Widaman, Helm, Castro-Schilo, Pluess, Stallings and Belsky2012) contrasted strong and weak versions of differential susceptibility theory predictions and advocated comparative testing of these versions. Under the strong version, the Environment main effect should be zero in persons with no plasticity alleles; under the weak version, the Environment main effect might be non-zero even in persons with no plasticity alleles. To test whether the strong version was supported in this data set, I restricted the Environment main effect to zero in Model 3. As seen in Table 7, the drop in R ² was miniscule, ΔR ² = −.001, and not significant, F(1, 276) = 0.33, p = .56, and the regression coefficient for the G×E interaction, B = 0.07 (SE = 0.01), was significant, t(276) = 5.45, p < .0001. The improved level of statistical significance for the G×E interaction was the result of the much smaller SE for the interaction coefficient in Model 3 (SE = 0.01) relative to the comparable value in Model 2 (SE = 0.05). As seen in Table 7, the test of the restriction of the E main effect to zero did not approach significance, t(276) = − 0.58.

I fit one final model, Model 4, which added a restriction that the G main effect was zero. Because the E main effect was mean-centered, this represents a test that the G effect was zero at the mean of the E predictor. As seen in Table 7, the drop in explained variance moving from Model 3 to Model 4 was small, ΔR ² = −.0006, and not significant, F(1, 277) = 0.23, p = .63, and the test of the restriction was nonsignificant, t(277) = − 0.48. The regression coefficient for the G×E interaction remained unchanged in Model 4, a model in which the G×E interaction explained about half of the 17 percent of variance explained by the equation. The adjusted R ² value for Model 4 was higher than comparable values for other models, and the BIC value was lower than for other models, attesting to the efficiency and optimality of Model 4.

Interrogating the weighting of SNPs

Models 1 through 4 utilized an equally weighted sum of the 5 SNP scores, so implicitly assumed that SNPs did not differ in their effects. This approach may have masked differential effects of the SNPs, so four additional regression models were tested, using the 5 dominant and 5 recessive SNP scores as separate predictors. In these analyses, the Gene main effect was represented by 10 separate dichotomous predictors, and the G×E interaction was represented by 10 product vectors, where each product vector was the product of a SNP score and parent hostility. With this approach, it becomes possible to interrogate or test more severely the equal weighting of SNPs.

The first of these models, identified as Model 5 in Table 8, used Female as control variable, parent hostility as Environment variable, and the 10 SNP scores to represent the main effect of Genes; this model explained over 17 percent of the variance in the outcome variable, R ² = .1755. As in the comparable Model 1 (Table 7), the E main effect was significant, B = 0.31 (SE = 0.07), p < .0001, but none of the 10 SNP scores met the .05 level of significance.

Table 8. Alternative Regression Models for Hostility using Individual SNP Codes for Masarik et al. (Reference Masarik, Conger, Donnellan, Stallings, Martin, Schofield, Neppl, Scaramella, Smolen and Widaman2014) Data

Note: N = 281. ^a The t ratio has df equal to ν ₂ for the F ratio for the equation, shown below. ^b For the F ratio for each model, ν ₁ = degrees of freedom numerator, and ν ₂ = degrees of freedom error (or denominator).

Model 6 added the 10 product vectors comprising the G×E interaction (see Table 8). Compared with Model 5, Model 6 explained more variance, R ² = .2128, but this increase in explained variance, ΔR ² = .0373, was nonsignificant, F(10,258) = 1.22, p = .28. Inspection of parameter estimates and their SEs in Table 8 reveals that not one of the G main effects had an effect significant at p < .05, and only one of the 10 G×E product terms (for the dominant effect of DAT) barely met the .05 level. In short, Model 6 looks like a veritable sea of nonsignificance.

However, restrictions on model parameter estimates can test whether coefficients differed from values that were predicted a priori. Model 7 invoked three restrictions: (a) due to presence of estimates for a G×E interaction, the main effect of the E variable was fixed at zero, (b) coefficients for the 10 gene SNP main effects were constrained to equality, and (c) coefficients for the 10 G×E interaction vectors were constrained to equality. In Table 8, Model 7, with R ² = .1724, explained less variance than Model 6, but the drop in fit, ΔR ² = −.0404, was not significant, F(19, 258) = 0.70, p = .82. This F-ratio is an overall or setwise test of the 19 restrictions in Model 7 and supports the contention that the 19 restrictions did no harm overall (or as a set). To ensure that no harm was done in particular, tests of parameter restrictions in Model 7 are shown in Table 9. There, the first restriction is the restriction of the E main effect to zero, Restrictions 2 – 10 reflect the equality constraint on the 10 SNP main effect coefficients, and Restrictions 11 – 19 reflect the equality constraint on the 10 G×E interaction coefficients. Not one of the tests of individual restrictions was significant at the .05 level

Table 9. Tests Restrictions on Regression Parameters in Models 7 and 8 for Masarik et al. (Reference Masarik, Conger, Donnellan, Stallings, Martin, Schofield, Neppl, Scaramella, Smolen and Widaman2014) Data

Note: N = 281. The beta distribution was used to compute probability levels for all t ratios (see SAS 9.4). ^a The t ratios in this column have 277 df; all ps > .15. ^b The t ratios in this column have 278 df; all ps > .17.

The final model, Model 8, added to Model 7 a constraint that the gene main effect coefficients be fixed at zero (see Table 8). Model 8 had a slightly lower level of explained variance, a drop in explained variance, ΔR ² = −.0006, that was neither large nor significant, F(1, 277) = 0.23, p = .63. To ensure that no harm was done to any single estimate, tests of restrictions in Model 8 are presented in Table 9. Not one of the 20 restrictions in Model 8 neared statistical significance, all ps > .22, so no single estimate was unduly harmed.

Note that Models 7 and 8 in Table 8 are identical in all ways – explained variance, parameter estimates and their SEs, adjusted R ² – to Models 3 and 4, respectively (cf. Table 7). Models 3 and 4 employed the Gene total score, so implicitly forced equality constraints on main and interactive effects of the gene SNP variables assuming that gene effects were (a) linear for each SNP and (b) equal in magnitude across SNPs. Models 7 and 8 enforced these constraints explicitly. The mathematical identity of all characteristics of Models 3 and 7 and of Models 4 and 8 is a brute force demonstration that the same end is achieved if one employs an a priori equally weighted sum of a set of indicators versus if one employs the indicators as separate predictors but forces equality constraints that embody those used in the a priori summation. Model 7, with its 19 restrictions, arrived at the same point as Model 3, but allowed the data to “complain” if any of the individual restrictions harmed the ability of any one of the genetic main or interaction components to contribute to the model. In effect, tests of restrictions in Model 7 allowed individual SNP components to “call foul” if a given constraint was overly restrictive. That not a single restriction was associated with statistical significance is a telling outcome that justifies the restrictions on a one-by-one basis. Similar comments apply to the identity between Models 4 and 8, as both of these models added the constraint that the gene main effect was zero.

In sum, Model 8 satisfies the three principles laid out early in the paper. First, Model 8 did no harm in general. It explained less variance than the most highly parameterized model, Model 6, but decreases in explained variance were nonsignificant. Second, Model 8 did no harm in particular, as not one of the 20 restrictions imposed on Model 8 led to any indication that a restriction harmed the contribution of any individual model component. Third, Model 8 did considerable good, with much smaller SEs of estimates relative to Model 6, a larger adjusted R ² than Models 5, 6, and 7, the lowest (i.e., best) BIC value of any of the models, and explicit justification of implicit assumptions made by Model 4. Thus, Model 8 reflects optimal return on investment when estimating regression parameters (cf. Rodgers, Reference Rodgers2019).

Structural equation modeling

I conducted comparable analyses to those reported in Table 8 using SEM programs Mplus (Muthén & Muthén, Reference Muthén and Muthén1998-2019) and the lavaan package (Rosseel, Reference Rosseel2012) in R. Essentially identical results were obtained; given limitations of space, the full set of analysis scripts and descriptions of SEM results were placed in Supplementary Material.