Introduction

1At the request of Daniel Courgeau and later Isabelle Séguy and Luc Buchet, I was drawn to a statistical question raised in paleodemography. I discovered a complex field of research, based on observations that are often few in number and of dubious reliability, often providing only indirect information on the issues studied: a challenge for the statistician (although far from the “big data” that is more fashionable today). I also discovered that some paleodemographers had highlighted the many challenges and paved the way, tackling some notorious errors, insisting on refraining from exaggerated enthusiasm, while looking for tools that would justify a more rational optimism. Jean-Pierre Bocquet-Appel is among them and has a prominent place in reflecting on the estimation of the adult age at death distribution based on aging biological data. The greatest concern for rigor is always present in Jean-Pierre’s work, even if the results lead to some disappointments. In statistics, rigor means a reliable method (without methodological bias) and an indication of the estimates’ precision level, considering the random nature of the observations. In addition, rigor requires maximizing efficiency by reducing this uncertainty inherent to all estimates as much as possible, without neglecting to broaden the scope of the method application. We will discuss Jean-Pierre Bocquet-Appel’s contribution to the implementation of such a program. We will present some recent works, mainly extending his reflections, and suggest how to pursue it further. It should be emphasized that our ambition is not to make a complete survey of the literature or compare in detail the different schools’ perspectives at any given stage of research development (see, for example, Hoppa and Vaupel 2002; Courgeau 2013). Nevertheless, we hope to cover the essential points around Jean-Pierre Bocquet-Appel’s fundamental reflection and technical contributions.

2We will arbitrarily begin this history at the beginning of the 1980s, with the following citation from Bocquet-Appel and Masset 1982 clearly outlining the preceding context:

“Historical demography […] succeeded in constructing a picture of life and death among the European peoples of three centuries ago […]. The contrast between the mortality pattern thus obtained and the one arising from our graveyard populations (some of them barely older) was striking. It challenged the quality of the conclusions propounded by the paleodemographers: were the populations they had been investigating actually as peculiar as they seemed to be? Or perhaps were the methods in use warped by some systematic error? The second hypothesis proved to be the right one. In fact, there are several causes of error, of a statistical rather than a biological nature. These causes of error were described in previous works (Masset 1971) […]. Let us only record the most awkward of them: the procedure which leads one to superimpose on the mortality structure of cemetery populations the structure of other populations entirely alien to them.”

3There are two main points to be emphasized here: the French school was quick to denounce the severe errors that marred several studies, and the statistical method must be challenged primarily. In order to understand the issue, it is essential to state the problem clearly, i.e., the question raised and the “raw material” available to try to answer it. The latter is of two kinds: available observational data and the halo of prior knowledge on the subject. The question under consideration is the distribution of ages at death on a site, for which it will therefore be necessary to forge the tools to move reliably and efficiently from the available means to the expected answer.

Age estimation from biological data: formalization of the statistical problem

The data

4We will present the problem for discrete (or discretized) data, i.e., defined by classes pertinent and often warranted in Jean-Pierre Bocquet-Appel’s work. Most of the biological aging indicators are particularly difficult to evaluate in a very precise way, and we know how much chronological ages evade a valid estimation to the extent that it would be foolhardy to expect an estimation within a few months. Moreover, the case where several aging biological indicators are simultaneously available will remain outside our scope. The required statistical methods are undoubtedly different in nature, in particular, because the valuable reference data are not easy to acquire (unless we assume the indicators’ independence, which is certainly inappropriate in most cases).

5Let r be the number of classes (stages) used for the indicator and c the number of age classes. The available data are:

• the frequencies of the indicator’s various stages for a given site, i.e., for a sample drawn from the population corresponding to that site; these frequencies will be noted: m₁, m₂, …, m_i,...,..., m_r, the sum of which m is the observed sample size,
• the frequencies observed in a reference table drawn from a population in which both the indicator stage and the age at death are known for each subject. We will write n_ij the number of subjects in the i^th skeletal stage and the j^th age class. Here again, we refer to frequencies observed in a sample (or rather in c samples, noting that it is always possible to compose a sample for each age class without modifying the various implemented methods; this is inherent in the methods themselves).

6Only the statistical issue will be addressed here. Therefore, we will assume that the important preliminary question of the quality of the observations has been resolved. For example, the guarantee (to a reasonable degree) that the collected samples represent the corresponding populations, that the coding of the biological indicator is analogous for the reference data and the data collected at the site of study, among others. However, one point should be mentioned regarding the reference data. The data are presumed to have the same distribution of skeletal stages for a given age in the reference population and in the site studied (the hypothesis of biological invariance implies that the reference data belong to the preindustrial period). On the other hand, for the age distribution of a given skeletal stage to be the same in the reference population and the site under consideration, it would be necessary to have an equivalent of biological invariance by permuting the roles of skeletal stage and age, which obviously cannot be achieved. This is precisely because the age structure of the site under investigation is unknown and has nothing to do in theory with the age structure of the reference population. This is another way of expressing the essential point of Bocquet-Appel and Masset’s quote from 1982 mentioned in the introduction, where the authors denounce “the procedure which leads one to superimpose on the mortality structure of cemetery populations the structure of other populations entirely alien to them.”

7In addition to the statistical data mentioned before, it is necessary to add the knowledge of the subject previously acquired by the researchers. All this prior knowledge can (and must) improve the statistical study.

Objective

8The goals are:

a) the estimation of the age structure of the (adult) population corresponding to the burial site studied, i.e., the estimation of the probabilities p₁, p₂, ..., p_j, ..., p_c of the c age classes in this population. However, the estimation of other parameters should also be considered, for example, estimating the average age at death for the adult individuals of this population.

b) the estimation of the level of uncertainty of the previous estimators.

c) the search for a statistical method (as general as possible) that reduces this uncertainty as much as possible.

9Before discussing how the methods have evolved to achieve the three objectives above successively, some preliminary notations should be introduced.

10We will denote π₁, π₂, …, π_i, ..., ..., π_r the probabilities of the various skeletal stages for the site under consideration. It is conventional for these parameters to be estimated by m₁/m, m₂/m, …, m_i/m, ..., m_r/m, respectively. We will note p_1/j, p_2/j, …, p_i/j..., ..., p_r/j the probabilities that the skeletal stage is respectively stage 1, stage 2, ..., stage r, for a subject in age class j. These probabilities are naturally estimated by n_1j /n_j, n_2j /n_j,...,n_rj /n_j, where is the sum of frequencies n_ij over all stages i (sample size of subjects of age j in the reference population).

From point estimation to uncertainty assessment

11Bocquet-Appel and Masset 1996 emphasize the importance of going beyond the search for a simple point estimator (valid on average) of the age structure by providing an assessment of its precision. To do so, they achieved a major breakthrough by proposing the use of the bootstrap technique, which had only been introduced into the arsenal of statistical tools recently (Efron 1979). This is a resampling method that appropriately varies the reference data around the observed values to simulate sampling errors on these data and hence on the estimation results. Bocquet-Appel and Masset compared the accuracy of their method, introduced long ago by C. Masset (IPFP algorithm) for point estimation, and the IALK algorithm proposed by Konigsberg and Frankenberg 1992 (following the work of Kimura and Chikuni 1987) on a fisheries problem formally similar to the one faced by paleodemographers. The competing algorithms IPFP and IALK have ultimately been proven identical (Bocquet-Appel and Bacro 1997). The minor differences sometimes observed are only approximations related to the practical implementation of the calculations1.

12Several years later, Konigsberg and Frankenberg 2002 insisted on the need for precision analysis, noting that the IALK algorithm provided a maximum likelihood estimate and proposed using the general properties of such estimators. However, these properties are only valid for sufficiently large samples, which is rarely the case in practice (for small samples, evaluating the precision by bootstrap would be more reliable).

13However, we notice that, apart from this sample size issue, either previous approach is still insufficient for a valid evaluation of the uncertainty. The data are effectively only samples allowing to approximate respectively the parameters π₁, π₂, …, π_i, ..., ..., π_r (site data) and, for each j, the parameters p_1/j, p_2/j, …, p_i/j..., ..., p_r/j (reference data). Therefore, there are two sources of uncertainty in any estimate using these data. However, Bocquet-Appel and Masset 1996 “bootstrapped” only the reference data without considering the variability of the site frequencies. Instead, the IALK method assumes that the reference data estimate the corresponding conditional probabilities without error, only the uncertainty of the site sampling being considered.

Extension of estimation possibilities

14The algorithm discussed above requires that the number r of stages is greater than or equal to the number c of age classes. This condition is easily understood with the IALK format based on a statistical model: if r<c the number of estimated parameters is greater than the number of available data, then the model is said to be unidentifiable (which is sometimes translated in the literature to mean that the number of degrees of freedom is negative!).

15One way to reduce the number of parameters is to introduce a demographic model for the p_j, thereby reducing the number of parameters from c-1 to the number of parameters of the model (e.g. 2 for the Gompertz model, 3 for the Gompertz-Makeham model). Specifically, with the Gompertz model and the minimum age of a years, the probability of the age class [b,c] years will be

depending on the two parameters λ and ρ alone regardless of the number of age classes. Several papers mention this possibility in the referenced volume by Hoppa and Vaupel 2002. Konigsberg and Herrmann 2002 proceeded to implement it. We include this option in Appendix 3 as part of a brief comparative study.

16With a similar objective, Jean-Pierre Bocquet-Appel introduces an innovative approach at this stage (Bocquet-Appel and Bacro 2008). Rather than using a model as in classical statistics, the authors seek an estimate in a family of mortality laws that they produce from prior knowledge derived from demographic data. The authors refer to this method as iterage, which can be described in the following way.

17The parameters of the model are connected by the equation:

(1)

18If p_i/j and π_i are replaced by their natural estimates, n_i/j /n_j and m_i /m, respectively, provided by the reference and site data of interest, the system of equations (1) becomes a multiple linear regression from which to estimate the parameters p_j. However, for c>r, there are too many unknown parameters, and the system has an infinity of solutions resulting in an ineffective method (this is the problem already mentioned above). Even if c<r, the usual least-squares estimation often gives odd or even aberrant results such as a negative probability estimate. Obviously, we can proceed by least squares under the constraints necessary to produce an estimate in the range of possible values for the probabilities p_j; however, we can also produce estimates on the range boundary, effectively a zero probability for some age groups, which is clearly unsatisfactory. This statistical method could probably be further improved on a theoretical basis. However, the subtle answer is to say: “the p_j cannot be arbitrary since it is a mortality law; consequently, we will seek the mortality law which best solves the system (1) (for example, using the least-squares method) among a family of potential candidate laws.” Moreover, since the reference data are subject to uncertainty, they can be bootstrapped around the observed values, resulting in as many estimated values as bootstrap repetitions. The mean of these values will provide a point estimate. In contrast, their distribution around this value will evaluate the uncertainty, measured by a standard deviation or any other dispersion index. The authors propose two families of “candidate” mortality laws, one for attrition and one for catastrophe.

19This method is effective in many circumstances, but two points need to be closely revisited:

• The two suggested mortality law families seem natural, but how to choose between them in a specific case? Should we sometimes combine them? Should we build others?
• For uncertainty calculation it is insufficient to bootstrap the reference data when the site data is itself a sample. Thus, the π_i estimate is also subject to uncertainty. Besides, the assessment of uncertainty presents a slightly more delicate challenge: the family of mortality laws considered is bounded, i.e., “corseted” between “extreme” laws; if the (unknown) law we are looking for is outside this family (more accurately outside the laws that can be derived by combining the members of this family), the solution found for each bootstrap replication will be nearest to (but possibly “distant” from) the same extreme value of the family or some of these extreme values, so that the dispersion of these solutions will be small. Therefore, the point estimate (average of these solutions) will be necessarily inaccurate with the illusion of minimal uncertainty (see Appendix 3).

More recent work

20The iterage method emphasizes the value of invoking general knowledge about the nature of the estimated parameters, particularly that p_j represents a mortality law. In practice, some parameters are likely to be observed more often than others. For a statistician, this implies the use of a Bayesian method. The methods of classical statistics (or frequentist) consider that the unknown parameter (in fact a vector of several real parameters, in our case, the p_j and the p_i/j) is fixed in a given space (the parametric space) and that the estimate is a value of this space supposed to approximate it. Specifying a model implies reducing the dimension of the parametric space (for example, by taking a Gompertz mortality model in the framework of the IALK method). In contrast, Bayesian statistics assigns a probability distribution to this parameter on the parametric space to weigh the different possibilities considering the general knowledge on the subject and independently of the observed data. This is why this distribution is called a prior law. The observed data indicate the specificity, on the one hand, of the joint distribution of the skeletal index and the age (reference data), on the other hand, of the site of interest (site data). They prompt to revise the vector distribution of the parameters accordingly: this distribution, conditional to the data, is called the posterior law of the parameters of interest2.

21In an attempt to improve the previous methods, we have proposed and implemented a Bayesian estimation method. An initial version of the method is outlined in Caussinus and Courgeau 2010 and Séguy and Buchet (chapter 13, 2013). This method is further developed and applied to new examples in Séguy et al. 2013. This article demonstrates that the number of age classes can be relatively high without damage (on the contrary): we recommend 13 age classes (one class of 18 to 24 years old and one class of over 80 framing 11 five-year classes). An improvement regarding the prior distribution of the p_j is outlined in Caussinus et al. 2017 and applied to new data in Buchet et al. 2017.

22We refer to the articles mentioned above for details of the method’s implementation. We will simply emphasize below the two fundamental aspects: the choice of prior laws and the nature of the results generated.

Prior Laws

23With the reference data, the samples we are working on are the c samples that provide each age class the frequency of skeletal stages. Their distribution is assumed to be multinomial with parameters p_1/j, ..., p_i/j, ..., p_r/j for the j^th one. There is no information on these parameters other than provided by the reference data. Therefore, the prior law must be “neutral”, which implies adopting a uniform law, a special case of the so-called Dirichlet laws (see Séguy and Buchet 2013, chapter 13, box 13.1) where all parameters are equal to 1.

24The situation is more complex for the prior distribution of the parameters p₁, …, p_j, ..., p_c. A first approach can be to transpose the iterage principle in Bayesian terms and to consider a prior uniform distribution on a family of mortality laws which could be the one constructed by Bocquet-Appel and Bacro 2008 or the “real” laws compiled by Séguy and Buchet 2013. This choice was not the most convincing one (see the referenced articles) but would certainly deserve to be considered again.

25A second approach is to start from a more theoretical perspective which leads to select a Dirichlet law since we are dealing (as above) with c positive numbers summing to 1. Since the p_j represent one mortality law among all possible laws, the selected Dirichlet law will be centered on a mean mortality law. To do so, we have opted for the preindustrial standard (Séguy and Buchet 2013).

26Subsequently, we sought to account for the correlation characteristics specific to all mortality laws. Irrespective of how these are considered, the variations in p_j from one class to another are relatively “smooth”, which led to the choice of a “smoothed Dirichlet distribution” as the prior law (constructed by a moving average on a Dirichlet distribution, always keeping the preindustrial standard as the mean: see Caussinus et al. 2017. This latter prior law will be used in the following examples.

Nature of the results produced - examples

27The calculations are performed by Monte-Carlo simulation. The desired estimates are derived from the posterior probability law of the parameters of interest. The mean (or median) of this law can be considered as a point estimate (central), its density also providing the degree of uncertainty; this can be quantified by its standard deviation or (rather) by quantiles: for example, the interval between the quantiles 0.025 and 0.975 of the posterior law is an interval, known as the 95% credibility interval, which plays the role of the 95% confidence interval in classical statistics. This procedure can be implemented for all parameters of interest, such as the various p_j, which is the first goal of the analysis, and for any function of these parameters. For example:

• the survival probabilities for s= 2, ..., c
• the mean age at death (of a subject who has reached adulthood) (where a_j is the center of the j^th age class), an approximate calculation that requires classes of sufficiently small amplitude,
• individual conditional probabilities (probability that an individual of stage i is in age class j).

28The following figures show some results with various representational methods on the Frénouville, Merovingian period data (for a presentation of this site data and the reference data used, see Buchet et al. 2017). The evaluation of aging is done using the cranial suture synostosis adjusted to five stages. Thirteen age classes are considered, namely 11 five-year classes between 25 and 80 years of age framed by one 18-24-year-old class and one over 80 (see Séguy et al. 2013 for a justification of this division).

1. Posterior densities of parameters p₁ (black) and p₁₂ (red) for the Frénouville site, Merovingian period.

29Figure 1 shows the posterior densities of the parameters p₁ (probability of the 18-24 age group) and p₁₂ (probability of the 75-79 age group). The central values (point estimates by the posterior average) are respectively 0.27 and 0.06 (indicated by the central vertical line for each density); the tighter pattern of the density suggests that the precision is significantly higher for the second probability than for the first. The quantiles of the posterior laws are a useful summary: the lateral lines specify the 2.5% and 97.5% quantiles, and therefore the interval between them is a 95% credibility interval.

2. Frénouville site, point estimate (+) of the 13 age class probabilities for the Merovingian and Gallo-Roman periods (shifted to the right), with credibility intervals at 50% (red segments) and 95% (green segments). Stars indicate the preindustrial standard probabilities.

30Figure 2 shows the point estimates (+ sign) for each age class with 95% (green between the 2.5% and 97.5% quantiles) and 50% (red between the 25% and 75% quantiles) credibility intervals. The graph provides a synthetic view of the potential differences between the Gallo-Roman and Merovingian periods and compares the preindustrial standard. The Merovingian period shows a clear difference with this standard (over-mortality of the first three classes), while the Gallo-Roman period is relatively similar. However, it is difficult to draw a firm conclusion on a potential difference between the two periods because of the overlap in credibility intervals. This difference is clearer when the age classes are grouped (Figure 3) and will be quite clear through the mean age at death.

3. Frénouville site, point estimate (+) of the three age classes 18-39 years, 40-59 years, and over or equal to 60 years for the Merovingian and Gallo-Roman periods (shifted to the right), with credibility intervals at 50% (red segments) and 95% (green segments)

31Figure 4 illustrates the posterior densities of the average age at death for Merovingian Frénouville (black) and Gallo-Roman Frénouville (red). However, calculating a mean age from an age class distribution is approximative, acceptable only for a sufficiently fine subdivision, as seems to be the case with the five-year classes. Moreover, the approximation due to discretization is the same for the two sites and will not affect their comparison. The tighter density for the Merovingian period corresponds to a larger sample size that allows for greater precision. Moreover, the two distributions displayed barely overlap, so that the average age at death is likely to be lower in the Merovingian period than in the Gallo-Roman period. The posterior probability of this event can be estimated from the densities by period and from the two site samples’ independence: the probability is estimated at 0.997, which is very close to 1, leading to the conclusion that the two mortality laws are different with great confidence.

4. Frénouville site, posterior probability densities of the average age at death for the Merovingian period (black) and the Gallo-Roman period (red).

32Other parameter functions are considered in Appendices 1 and 2.

Conclusion and perspectives

33About thirty years ago, Claude Masset and Jean-Pierre Bocquet-Appel finally convinced paleodemographers that their perspectives on the age estimation based on an aging indicator were valid.

34The fundamental developments on this question, such as evaluating the estimates’ precision, expanding the methods to broader situations, and improving their efficiency, could begin. Jean-Pierre Bocquet-Appel was present at each step of this process. With the iterage algorithm, Jean-Pierre Bocquet-Appel provided the transition between the frequentist statistical methods exclusively used so far and the Bayesian methods promoted since. We have discussed (see appendix 3, in particular) some challenges of the iterage algorithm, such as the choice of “candidate” set mortality laws and boundary problems linked to the set necessarily abrupt edges.

35The prior distribution introduced in the Bayesian methods does not present these edge problems, but does it sufficiently consider all the aspects of a mortality law? In competition with iterage, Bayesian methods sometimes seem to have lower accuracy. Is this detrimental or simply realistic? In other words, is the higher accuracy illusory or effective? Clearly, both cases are possible, but we never know in which situation we operate. This suggests that the Bayesian methods’ prior law could be improved in light of Jean-Pierre Bocquet-Appel’s proposals. Although some attempts have been made by constructing a prior law on a family of mortality laws in a simple way (see, for example, Caussinus and Courgeau 2010, 2013; Caussinus et al. 2017), all still have drawbacks, and other avenues should be explored. For instance, the proposed Dirichlet law smoothing also follows the same logic but remains very heuristic. How could it be better systematized? How could it benefit from considering lists of mortality laws?

Appendix 1. Note on the results’ interpretation from a Bayesian analysis

A1. Point estimates of posterior survival probabilities for the Frénouville site, Merovingian period: males (M), females (F), and sexes combined (T). Stars indicate survival probabilities of the preindustrial standard for all sexes.

36The principle of Bayesian statistics is to investigate how the observations require a revision of the prior distribution. Therefore, this prior law plays a pivotal role in interpreting the results of any statistical analysis. We can illustrate this with an example. Figure A1 shows the point estimate of survival probabilities for men, women, and both sexes combined at the Frénouville site (Merovingian period). The survival curve corresponding to the preindustrial standard has been added. We might be surprised that the survival curve for both sexes combined is not just an average of the survivals by sex. However, we must remember that these estimates are always derived by adjusting the prior law according to the observations. As the prior law is centered on the preindustrial standard, if the population age structure is far from the standard, the prior law will be adjusted. However, this adjustment will be increasingly important as the data become more informative, i.e., more abundant. Moreover, the sex-specific sample sizes (at best about half of the total, in this case, 61 and 74 against 200) result in a significantly lower level of information than the total sample, which explains the paradox3. Additionally, the figure indicates a lower estimated survival for women, but the sex-specific credibility intervals (not shown) overlap substantially, suggesting that the apparent difference is not necessarily significant. A synthetic parameter such as the mean age at death could be more discriminating. We find a point estimate of 47.85 years for men and 43.72 years for women. The posterior probability that this age is higher for men is 0.812. This value is quite high but still below the levels usually required for a definite conclusion, especially if we use the traditional “level” of the significance test paradigm in the present framework. Nevertheless, from an elementary classical statistic perspective, we observe that a test of homogeneity (chi-square) between the frequencies of the stages in men and women indicates a probability of excess (p-value) of 0.09, not significant at the usual levels. However, the test suggests a strong tendency towards a marked difference between the sexes (although without necessarily explaining why). Somehow, this echoes a point of view that reflects Jean Pierre Bocquet-Appel’s (and Claude Masset’s) thinking: in paleodemography, the nature and parsimony of the information that can be gathered must lead to the greatest caution and, while it is advisable to seek out the most suitable methods, one must not “believe in miracles”.

Appendix 2. Individual age

37Although it may have been tempting to try estimating the age at death of each individual at a site from its skeletal stage and properly aggregate these results to estimate the age structure of the site, this approach was long ago proven to be ineffective. Paleodemographers have therefore considered estimating the age structure directly. However, as mentioned previously, the probability that an individual of stage i is of age class j is a function of the parameters p_i/j and p_j that can be estimated within the framework of the procedure in paragraph 5. This is easily extended to estimate the probability that an individual of stage i is of an age falling into any of several classes. To illustrate this, Figure A2 displays the posterior densities of four of these parameters for the Frénouville site (Merovingian period), and Table A2 gives the point estimates and the 95% credibility intervals. First, we observe that the probability of a single age group can be relatively small, even if it is an “extreme” age associated with an “extreme” stage. The probability is obviously higher for a larger age group, although, for a fairly large “central” age group and a “central” stage, the probability turns out to be hopelessly small. Most importantly, the 95% credibility intervals are wide (their length is between 0.228 and 0.355), indicating very imprecise estimates. This result reflects the difficulty of estimating this type of conditional probability correctly. Nevertheless, if it is helpful to estimate conditional probabilities in some cases, it can be done but with great precaution... and the risk of disappointing results.

A2. Frénouville site, Merovingian period, posterior density of some conditional probabilities with .025 and .975 quantiles.

Appendix 3. Comparative elements of various approaches

38We have encountered fundamentally different statistical methods in reviewing the methods’ evolution for estimating an adult age structure. At the outset, we found the equivalent approaches of ITFP and “basic” IALK. The latter can be extended to any number of age classes by considering a demographic model. The method is based on classical frequentist statistics. The model allows to reduce the number of parameters, i.e., the dimension of the parametric space, and the estimation is done using maximum likelihood. Instead, we have presented a Bayesian approach in which the probabilization of the parametric space replaces the drastic reduction of the frequentist method. The iterage method of Bocquet-Appel and Bacro 2008 is intermediate (in time and essence) between the two previous ones. We will compare these three methods on an example where register data exist. The Antibes site (19th c.) is presented in Séguy et al. 2013. The biological indicator is the cranial suture synostosis, for which five stages have been selected. We divided the ages into eleven five-year classes, framed by the first class from 18 to 24 years of age and the last class over 80. Since the prior distributions provided by Jean-Pierre Bocquet-Appel were not adapted to this age distribution, we opted instead for the 167 mortality tables compiled by Séguy and Buchet 2013 (we can check their resemblance to the “attri” distributions of Bocquet-Appel). Furthermore, we fitted a Gompertz model on the results provided by methods that do not explicitly use it (nonlinear least-squares fitting). Figure A3.1 and Table A3 reveal that iterage and the Bayesian method perform very closely and significantly better than the frequentist method, including when searching for the model that the latter implies. The point here is not to set up an efficiency contest on a singular test that may be biased (for that, a large comparison would have to be organized, probably by simulation, which is outside our scope). However, the consequent results raise interesting questions that deserve to be explored in greater depth. A Gompertz model fits the register data quite well and still approximates the data well when fitted to Bayesian estimates. The same applies when the model is fitted to iterage estimates. So why doesn’t the maximum likelihood method based on this model provide a good estimate? Is this an intrinsic benefit of the Bayesian method or a similar method, and if so, why? Could it be related to the similarity of the age structure of the site to the standard structure? (We checked that the choice of the model is not the issue: rather, the Gompertz model is the most appropriate here).

A3.1. Structure estimation using three statistical methods: Bayesian method from paragraph 5 (B - dark blue), iterage (I - red), IALK with Gompertz model (G - green). Also, we show the probabilities extracted from the registers (* - black), the Gompertz model fitted to the estimates by the Bayesian method (A - light blue), the Gompertz model fitted to the register data (R - pink).

39To examine some aspects of the comparative behavior of the iterage algorithm and the Bayesian method, we return to the example of Frénouville (Merovingian period). The age is divided into seven classes: five decennial classes from 30 to 80 years old framed by a “young adults” class and an “over 80” class, to use the two sets of “candidate” mortality laws provided by Jean-Pierre Bocquet-Appel: Attri for attritional and Cata for catastrophic. The iterage method is implemented with the modification suggested by our second remark at the end of paragraph 4: implementing a double bootstrap to account for sampling uncertainties on both the reference and the site data. The Bayesian method is applied with 13 classes which are then grouped into seven classes according to the present breakdown. Figure A3(2) shows the estimate of the age structure using the three methods discussed, each time with the estimated precision: 95% confidence interval (approximated by the rule of plus or minus two times the estimated standard deviation) for iterage, 95% credibility interval for the Bayesian method using quantile estimation. These intervals are globally longer with the Bayesian method, giving an intermediate point estimate between the two iterage estimates. The latter result might be expected to derive from the choice of candidate laws for iterage. However, if we apply iterage by concatenating the two sets of candidate laws, we arrive at roughly the same results as with the Cata set alone. Hence, the site seems to correspond to this situation (which incidentally raises the question of the prior choice between both Attri and Cata): perhaps it would be better to simply put them together instead. However, the estimation of the precision seems to raise the most interesting questions. Disregarding the fact that confidence intervals and credibility intervals are not quite equivalent, the overall shorter iterage intervals seem to argue for this method. However, by examining the set of candidate mortality laws, it appears that for solutions of system (1) with the various bootstrap replicates, the algorithm leads to the choice of laws close to the extremes, for example, in the 30-39 age group, a low extreme with Cata and a high extreme with Attri. However, we do not know which “true” mortality law is being estimated. If the law is at the center of the set of candidate laws, it will be accurately estimated, and the low dispersion warranted. If the law is outside this set, it will necessarily be poorly estimated. All the bootstrap replications will provide a candidate law as close as possible to this law outside their family. Thus, the laws will generally be close to each other, leading to a low variance and thus to the illusion of high precision.

A3.2. Estimation of the age structure for the Frénouville site, Merovingian period (circles), with estimated 95% confidence or credibility intervals (vertical segments): Bayesian method (black), iterage with Attri candidates (green), iterage with Cata candidates (red).