## Abstract

The identification of parity effects on the hazard of a next birth in cross-family data requires accounting for heterogeneity in fecundity across couples. In a previously published article, Cinnirella et al. *Demography, 54,* 413–436 (2017), we stratified duration models at the maternal level for this purpose and found that the hazard of a next birth decreases with rising parity in historical England. Clark and Cummins *Demography, 56* (2019) took issue with this finding, claiming that the result is a statistical artifact caused by stratification at the maternal level. This reply documents that our previous finding is robust to addressing Clark and Cummins’ critique.

## Introduction

Traditional scholarship has suggested that marital birth control was uncommon before the demographic transition. In a previous study (Cinnirella et al. 2017), we contested this view, showing that couples in pre-transitional England regulated their birth intervals in response to real-wage changes (parity-*independent* birth control) and to the number of surviving offspring (parity-*dependent* birth control). These results were based on duration analysis using individual-level statistics from the Cambridge Group’s *family reconstitution* data covering the period 1540–1850.

The comment by Clark and Cummins (2019), printed in this issue, suggests that our previous estimates of parity-*dependent* birth control are biased. Using simulated data, they claimed that the estimated coefficients of net parity are biased in models that stratify on the level of the mother. They also argued (without empirical substance) that we ignored important factors in our analysis of parity-*independent* birth control, claiming that the results could be “potentially explained as simply the product of nutrition, living environments, and social patterns.” Based on this reasoning, Clark and Cummins declared that “the traditional view of an absence of any evidence of birth control through spacing is vindicated.”

We are grateful to Clark and Cummins for replicating our earlier estimates and for giving us the opportunity to further examine the existence of birth control in pre-transitional England while accounting for their criticism. The analysis presented in this reply establishes that Clark and Cummins’ critique has no bearing on the conclusions presented in our earlier article (Cinnirella et al. 2017). In particular, the evidence of parity-*dependent* birth control is remarkably robust to the use of a range of identification strategies *not* subject to Clark and Cummins’ concerns. Although the magnitude of the estimates presented in this reply depends on the particular method used to account for the heterogeneity in fecundity across households, the conclusion that within-marriage birth control existed in pre-transitional England still holds. We also take the opportunity here to recapitulate the empirical strategies and key control variables we used in our earlier article to mitigate any potentially confounding factors in the analysis of parity-*independent* birth control. The present study thus solidifies the conclusion that pre-transitional couples engaged in parity-*dependent* as well as parity-*independent* birth control.

In particular, we estimate the effect of net parity on the hazard of a subsequent birth by using three model specifications that are not affected by Clark and Cummins’ concerns: (1) we estimate the relationship between *relative* net parity and birth spacing; (2) we account for between-family differences in fecundity, stratifying by the *protogenesic interval* (i.e., the period between the couple’s marriage and their first birth); and (3) we use duration models that account for family-level shared frailty. Even after the culprits that Clark and Cummins identified for our alleged statistical artifacts are eliminated, we still find statistically significant negative effects of net parity on the hazard of a next birth, further strengthening the evidence for the existence of pre-transitional marital birth control.

## Methodological Issues

As we explained in our earlier article (Cinnirella et al. 2017:415), more-fecund couples are more likely to be overrepresented at higher parities than less-fecund couples. Unobserved heterogeneity in couples’ fecundity (and other co-determinants of larger family size and shorter birth intervals) may thus generate an upward omitted-variable bias in the estimates of the association between parity and the hazard of a next birth. This upward bias may obscure parity-dependent fertility limitation in cross-family estimations when heterogeneity in couples’ fecundity is not accounted for.

We used Cox proportional hazard models, stratified on the level of the mother, in order to account for unobserved heterogeneity in couples’ fecundity (i.e., maternal-level fixed effects). Moreover, in order to control for the fact that maternal fertility is affected by age, we also controlled for the mother’s age at the beginning of each of her birth intervals. Clark and Cummins, however, asserted that the dual control for mother fixed effects and maternal age downwardly biases the estimates of net parity coefficients. Based on this assertion, they claimed that “[t]he parity effects are […] an artifact of the bizarre fit of the age controls to the data” (Clark and Cummins 2019). Importantly, they did not demonstrate their assertion in the context of a formal statistical framework. Rather, they loosely inferred their assertion from the outcome of a trivial simulation exercise. As a result, the size of the potential bias, as well as the conditions necessary for its existence, are not clearly identified. Nevertheless, the concern that they alluded to is potentially valid and thus deserves a proper consideration.

In the next section, we propose three alternative strategies, each of which accounts for the criticism Clark and Cummins raised and thus serves to examine the robustness of our earlier results (Cinnirella et al. 2017). Indeed, we document in this reply that net parity effects are also evident when we do *not* stratify on the level of the mother and account for heterogeneity of couple in the other ways. The sensitivity checks establish that the finding of birth control effects before the demographic transition is is robust to a wider range of additional statistical specifications than already established by Cinnirella et al. (2017).

## Alternative Models

In this section, we present further evidence corroborating the existence of net parity effects on the hazard of a subsequent birth in historical England. As discussed earlier, when identifying deliberate birth limitation through net parity effects, the researcher must account for heterogeneity in couples’ fecundity. Clark and Cummins (2019) claimed that mother-level stratification alongside the inclusion of control variables for maternal age at the beginning of the mother’s birth intervals prevents the production of “a reliable causal parity effect estimate.” Here, we propose three alternative models that account for heterogeneity in fecundity between couples but do not rely on mother-level stratification.

### Models With Relative Birth Order

The first specification assesses the correlations between birth spacing, parity, and family size. It demonstrates that net-parity coefficients smaller than one do not depend on the inclusion of mother fixed effects and controls for mother’s age at childbirth. In particular, we account for heterogeneity in couples’ fecundity by focusing on the *relative net parity*. By focusing on measures that capture the *relative* position of births among all family births, we provide a first attempt to disentangle family size from birth order, implicitly accounting for the unobserved heterogeneity in fecundity across families but *without* including mother fixed effects.

The relative net parity measure is defined as (*n* – 1) / (*N* – 1), where *n* is the net parity (i.e., the birth order of the surviving children) of the child that marks the beginning of the birth interval, and *N* is the total number of births in the family (see Ejrnæs and Pörtner 2004). Based on this measure, the first surviving child has a relative parity equal to 0, and the last surviving child has a relative net parity equal to 1. We also generate a series of dummy variables indicating the quartiles of the relative net parity measure in order to assess the impact of different values of relative net parity on birth intervals. The estimates of the duration models using the relative net parity are presented in Table 1. In column 1, we show the estimates produced by a continuous measure of relative net parity. The coefficient for relative net parity establishes that higher relative net parity is highly significantly associated with a lower hazard of a subsequent childbirth. In column 2, we include indicators for quartiles of the measure for relative net parity. The gradients of the coefficients document that the association between relative net parity and the spacing of births is monotonic: higher parities are related to a lower risk of having a next birth.

The estimated relationship between relative net parity and birth intervals could, however, be biased by the inclusion of the total household fertility in the relative birth order index because total fertility is potentially endogenous.^{1} Indeed, *ceteris paribus*, shorter birth intervals are expected to increase the final family size, resulting in a reverse-causality problem and thus potentially generating biased estimates. We attempt to shed light on this issue by controlling for net parity and family size separately. The results of this specification are presented in column 3. The coefficients for net parity are very similar to the relative coefficients reported in column 1, whereas total household fertility is associated with shorter birth intervals. Still, an endogeneity bias cannot be ruled out in the context of the present empirical strategy. Our empirical strategy in the models detailed further later in the article is therefore designed to avoid potential issues of endogeneity entirely.

Finally, as we discussed in detail in our earlier publication (Cinnirella et al. 2017:430–432), households drop out of the sample when mothers leave the parish in which they are observed. If birth spacing and sample attrition are systematically related, then this attrition could potentially bias our estimates. In particular, if longer birth intervals at higher parities are underrepresented because of sample attrition, then we might underestimate the existence of parity-dependent birth control. To avoid this issue, we constrain the sampled households to *completed* marriages—that is, those in which both spouses survive (in marriage) until the wife turns 50, thus ending her reproductive period. Estimates of relative net parity effects among *completed* marriages are reported in column 4. Indeed, the coefficients for relative net parity are lower than the estimates based on the full sample, suggesting that sample attrition might affect net parity coefficient estimates.

The concept of *relative* net parity helps us to disentangle parity progression from family size. The relative net parity specifications demonstrate that the observed net parity effects are due to neither the use of mother fixed effects nor some “impossible” mother’s age effect, as Clark and Cummins imagined. However, the use of relative net parity may generate an endogeneity bias of our coefficients of interest, as mentioned earlier in this reply. In what follows, we thus propose a set of alternative specifications that in various ways account for unobserved heterogeneity while avoiding issues of endogeneity.

### Models With the Protogenesic Interval

Our first alternative model, designed to estimate the effect of parity on birth spacing, makes use of the *protogenesic interval* as a proxy for a couple’s fecundity. The protogenesic interval captures the time elapsed between a couple’s marriage and their first birth. As Klemp and Weisdorf (2019) and Galor and Klemp (2019) detailed, longer protogenesic intervals are indicative of lower levels of parental fecundity among parents who did not conceive before marriage. Moreover, the protogenesic interval is highly significantly correlated with the length of subsequent birth intervals in historical populations (Galor and Klemp 2019; Klemp and Weisdorf 2019). Therefore, couples with shorter protogenesic intervals and hence higher levels of fecundity, *ceteris paribus*, are able to reach higher parities within a given time frame, as reflected by the fact that the protogenesic interval is a strong predictor of completed family size among historical populations. The negative association between the protogenesic interval and family size based on the Cambridge Group’s family reconstitution data is illustrated in Fig. 1.

Moreover, studies using modern data have suggested that the “time to pregnancy”—a concept closely related to the protogenesic interval—is not systematically related to unobserved potentially confounding factors (Klemp and Weisdorf 2019). Likewise, Klemp and Weisdorf (2019) documented, on the basis of the current historical data set, that the protogenesic interval is not statistically significantly associated with otherwise potentially confounding factors, such as parental socioeconomic characteristics (e.g., occupational skills and wealth).^{2} The protogenesic interval does contain a random component, which means that the interval is not a perfect proxy for fecundity. But controlling for the length of the protogenesic interval still mitigates any omitted variable bias associated with between-family heterogeneity in fecundity.

To this end, we restrict the sample to couples whose protogenesic intervals are at least 40 weeks long (i.e., the first child was conceived after marriage) and less than five years long (because long protogenesic intervals could be indicative of an unobserved first birth). Next, we group couples into deciles of the distribution of the protogenesic intervals and stratify our baseline Cox proportional hazard model on these groups. Importantly, we do *not*, as Clark and Cummins suggested, stratify on the mother level in this case. This specification then allows couples with different protogenesic intervals (and thus different levels of fecundity) to have different baseline hazards. In light of the plausible existence of an upward bias in the estimation of the effect of parity on the hazard of a birth in models that do not account for heterogeneity in couple fecundity, as described in the Methodological Issues section, our prior is that the estimated effect of net parity on the hazard of a next birth is consistently *lower* than when we control for couples’ protogenesic intervals.

The estimates of the duration models are detailed in Table 2. Column 1 reports the baseline coefficient estimate on the net parity variable when the model is not stratified by the protogenesic interval and therefore does not account for heterogeneity in fecundity. The specification stratifies by parish and quarter-century and includes a range of baseline control variables.^{3} The coefficient estimate for net parity is not significantly different from 1. In column 2, we show the results of estimating the same model but instead using dummy variables for net parity. The coefficients for net parity initially decrease but eventually increase with higher parities. Although the initial decrease is consistent with parity-dependent birth control, the subsequent increase suggests the existence of omitted variable bias. It concurs with the notion that the omitted variable bias associated with *not* accounting for between-family heterogeneity in fecundity can obscure a true underlying negative effect of parity on the subsequent hazard of a birth and is consistent with the findings presented in our earlier article (Cinnirella et al. 2017: Table 4, column 1).

One of the central concerns of Clark and Cummins is related to the estimated relationship between maternal age and birth spacing in the models that are stratified on the mother level. The coefficients for the quadratic polynomial in mother’s age in columns 1 and 2 of Table 2 suggest a convex relationship between age and the hazard of a birth. In Fig. 2, we report the predicted marginal effects of maternal age on the hazard ratio for the specification without stratification by protogenesic intervals (Models 1 and 2 in Fig. 2).^{4} The observed relationship is consistent with the notion that fecundity declines with maternal age. The inclusion of the control for the last birth interval generates a slight upward trend in the quadratic fit at the highest maternal ages, as revealed in Model 3 of Fig. 2, which displays the marginal effects of mother’s age at childbirth for a specification without a control for the last birth interval.

In column 3 of Table 2, we use a more flexible nonparametric specification, which includes dummy variables for mother’s age at birth in five-year intervals. The coefficients for the dummy variables lead to the same conclusion—namely, that birth hazards decrease with maternal age. As our new analysis establishes, there is nothing “impossible” about the estimated relationship between age and the hazard of a successive birth.

In columns 4 and 5 of Table 2, we show the estimates for models that are stratified by the protogenesic interval. A single net parity variable is used in column 4, and the net parity indicator variables are used in column 5. The estimates establish that the coefficient for net parity is significantly less than 1 after we account for differences in fecundity as proxied by the protogenesic interval. This reveals that higher parities are associated with longer birth intervals and thus lower fertility (column 4), a clear sign of parity-dependent fertility control. Likewise, when we use dummy variables for net parity, the coefficients tend to decrease with parity and stabilize at values below 1 (column 5).

Figure 3 reports the marginal effects of mother’s age at childbirth for different specifications for stratification by protogenesic intervals. The graphs report a total of four models: (1) the baseline specification with a continuous variable for net parity (Model 1; column 4 of Table 2); (2) the specification with net parity dummy variables (Model 2; column 5 of Table 2); (3) a modified version of the second specification that omits a dummy variable for the last birth interval (Model 3); and (4) a version that includes neither controls for net parity nor for the last interval (Model 4). Figure 3 shows that the hazard of a next birth declines with maternal age regardless of the specification. As for the specification without stratification on the protogenesic intervals, Model 3 shows that it is the inclusion of a dummy variable for the last birth interval that generates a slight upward slope in the quadratic fit toward the end of the reproductive period. Importantly, by excluding both controls for net parity and last birth interval, Model 4 establishes that the inclusion of the parity variable in the other models *does not* qualitatively affect the shape of the quadratic fit.

Because the use of a quadratic fit may be overly restrictive, we also exploit a more flexible specification with dummy variables for maternal age (column 6, Table 2). We find that the risk of a next birth decreases with net parity independently of the negative effect of age on birth risk, which strongly indicates the existence of parity-dependent birth control. Figure 4 reports the coefficients for the dummy variables of mother’s age at childbirth, estimated in columns 3 and 6 of Table 2. The maternal age–fecundity relationship estimated in models that stratify by the protogenesic interval does not present any “impossibility,” as shown in Figs. 3 and 4.

Finally, Fig. 5 provides a graphical representation of the coefficients for net parity reported in columns 3 and 6 of Table 2. The graph shows that the fecundity-adjusted coefficients for net parity are, as predicted, systematically smaller than the unadjusted coefficients and that they are systematically below 1. This set of results suggests that (1) accounting for fecundity mitigates the omitted variable bias associated with between-family heterogeneity, and (2) net parity has a negative effect on the hazard of a next birth. Again, these results for net parity are not obtained *because of* or *in the presence of* an “impossible” age-fecundity relationship, as Clark and Cummins claimed. These (new) results support the notion of a parity-dependent birth spacing effect and, by implication, are inconsistent with Clark and Cummins’ notion of random pre-transitional marital fertility.

### Models with Shared Frailty

Shared frailty models present yet another way to deal with Clark and Cummins’ concern. Duration models with shared frailty are effectively random-effects models in which, in their simplest form, an unobserved random factor modifies the hazard function of an individual or a group of individuals multiplicatively. Duration models with shared frailty are therefore useful to account for heterogeneity in couples’ fecundity *without* mother-level stratification. Models with shared frailty have been previously used in the historical demographic literature (e.g., Bengtsson and Dribe 2006; Dribe and Scalone 2010). A complication with the shared frailty model is its high computational cost: the model estimates random effects for each of the sampled households. Hence, we estimate shared frailty models on a 5 % and a 10 % randomly selected subsample of the data set.

The results of these specifications for both samples are reported in Table 3. In columns 1 and 3, we include a continuous variable for net parity; in columns 2 and 4, we include dummy variables for net parity. The relative coefficients for net parity from columns 2 and 4 are displayed in Fig. 6. Consistent with the previous results, we find that higher net parities have a monotonic negative effect on the hazards of subsequent births, further supporting the hypothesis of deliberate birth control. Also, as expected, maternal age has a positive association with the length of birth intervals. Importantly, we also find a significant frailty effect, suggesting the existence of potentially important unobserved heterogeneity across households.

## The Effect of the Real Wage on Birth Spacing

In a crude and highly speculative fashion, Clark and Cummins also attempted to raise doubts about the parity-*independent* analysis of birth control we reported in our earlier publication (Cinnirella et al. 2017). There, we presented evidence that couples responded to changes in real wages by regulating the length of their birth intervals. In particular, we showed that couples tended to postpone the birth of the successive child during periods of relatively low real wages, and vice versa. Clark and Cummins (2019) claimed that what we showed are just “connections,” that the relationship between occupational groups and the spacing of their births could “well be just a product of the location and lifestyles of that group,” and that the relationship between spacing and real wages “could be just a product of nutritional stress.” Given the arbitrary nature of these claims, along with a complete lack of any empirical demonstration that contradicts our findings, we suspect that their claims are based on a misreading of our article. In this section, we therefore briefly recapitulate how we accounted for potential confounding factors in our empirical analysis.

When studying the relationship between occupational groups and birth spacing, we stratify by parish and quarter-century. Stratification by parish accounts for time-invariant characteristics that vary across parishes, thus controlling for “location” or “environment” effects that Clark and Cummins mentioned. Stratification by quarter-century accounts for potential “lifestyles” that change with time and that affect all the parishes in the same way.^{5}

Moreover, we conducted an additional analysis on data restricted to six occupational groups while stratifying on the level of the mother and quarter-century (Cinnirella et al. 2017: Table 8). For all six groups, real wages were negatively correlated with birth intervals. When we combined the most affluent groups (i.e., farmers, merchants, and gentry), the real-wage coefficient was also statistically significant (Cinnirella et al. 2017: Table 8, columns 5 and 6; coefficient estimate = 1.073; *p* value = .049). These results again suggest that the negative association between real wages and birth intervals is *not* a product of differences in living environments or social patterns.

We also considered the effects among subsamples, including stayers, migrants, families with relatively short birth intervals, and *completed* marriages (Cinnirella et al. 2017: Table 7). The real-wage coefficient was significantly larger than 1 in all subsamples, suggesting that heterogeneity in environments and social patterns among these subgroups did not account for the observed *preventive-check* effect within marriages.

Finally, in an additional analysis (Cinnirella et al. 2017: online appendix, Table S6), we accounted for the potentially confounding effects of climatic conditions (proxied by temperatures) as well as nutrition and the living environments (proxied by crude death rates). The effect of real wages on birth intervals remained highly statistically significant. Overall, our earlier analysis (Cinnirella et al. 2017) clearly demonstrates that social patterns, living environments, and nutrition are unlikely to be accountable for the observed parity-independent birth control behavior.

## Conclusion

The identification of net parity effects on subsequent births using duration models requires accounting for heterogeneity in couples’ fecundity. In our previous article (Cinnirella et al. 2017), we used duration models that addressed heterogeneity by stratifying observations on the mother level. Based on family reconstitution data from historical England, our study revealed a negative and statistically significant effect of net parity on the risk of a subsequent birth. We interpreted this finding as evidence of birth control within marriage prior to the demographic transition. Clark and Cummins (2019) took issue with this finding, claiming that it is a statistical artifact caused by stratifying on the maternal level while controlling for the mother’s age at the beginning of each birth interval and concluding that marital fertility in pre-transitional England was random.

In this reply, we demonstrate that Clark and Cummins’ claim and conclusion are wrong. In particular, we show that the negative effect of net parity persists *even after* we account for heterogeneity in fecundity through stratification by the couples’ protogenesic intervals and the use of models with shared frailty. Although Clark and Cummins could have been correct in assuming that the joint mother-level stratification and inclusion of control variables for maternal age might have downwardly biased net parity coefficients, the existence of deliberate fertility control in the form of net parity effects remains evident. We show that estimates of the hazard of a next birth increases significantly with net parity independently of the identification strategy used. Moreover, Clark and Cummins were also wrong in asserting that our findings are dependent on “impossible” maternal age estimates. Maternal age is overall negatively related to the hazard of a next birth. Furthermore, in Clark and Cummins’ arbitrary criticism of our results showing parity-independent birth control, they ignored our extensive set of control variables, stratification strategies, and subsample analyses used to account for potential confounding factors.

Three important lessons spring from the additional analyses presented here. First, any serious attempt to identify deliberate birth control cannot ignore unobserved heterogeneity in fecundity across households. Second, the estimated magnitude of net parity effects on the risk of a next birth depends on the specific method used to account for such heterogeneity. Finally, regardless of the method used to account for unobserved heterogeneity, ample evidence of deliberate marital birth control exists for pre-demographic transition economies, including historical England.

## Notes

^{1}

The results presented in this subsection are robust to the use of a *predicted* measure of the final family size calculated for each birth interval by dividing the time remaining until the age at onset of age-related sterility by the average length of birth intervals. When the same constant values are used for age-related sterility and average length of birth interval for all families (e.g., sample averages), such a predicted measure is not affected by endogeneity concerns (related to, for example, an effect of birth spacing on final family size), which may potentially affect the observed measures.

^{2}

Furthermore, in our prior study (Cinnirella et al. 2017), we found almost identical effects of the real wage on the timing of the marriage date and the date of the first birth, indicating that the protogenesic interval (in contrast to interbirth intervals) is unrelated to the real wage rate and corroborating the notion that the marriage marked the intention to conceive. (Note that in our prior article, the text “a positive and statistically significant correlation between the real wage and the protogenesic interval” should instead read “a positive and statistically significant correlation between the real wage and the hazard of the first birth” (p. 420).) Indeed, a family-level Cox proportional hazard regression of the protogenesic interval on the real wage rate and the control variables (Cinnirella et al. 2017: Table 3, columns 1–2) yielded an estimate of the hazard rate of wages of 0.99 (*p* = .40).

^{3}

These baseline control variables include indicators for the occupation of the father in order to account for additional sources of relevant heterogeneity between families. Our results do not depend on the inclusion of these variables, however.

^{4}

The marginal effect of mother’s age is computed for the following: net parity 3; no mortality of the previous child; no last birth interval; occupation craftsmen; and not born on January 1, January 11, or December 25.

^{5}

Furthermore, the joint stratification by parish and quarter-century accounts for location or environment effects specific to each quarter-century.

## References

## Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.