{'fraction_non_alphanumeric': 0.0372341567053673, 'fraction_numerical': 0.009957722600167392, 'mean_word_length': 4.177555555555555, 'pattern_counts': {'":': 0, '<': 14, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 47, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Estimation of temperature-dependent growth profiles for the assessment of time of hatching in forensic entomology

Forensic entomology contributes important information to crime scene investigations. In this paper, we propose a method to estimate the hatching time of larvae (or maggots) based on their lengths, the temperature profile at the crime scene and experimental data on larval development. This requires the estimation of a time-dependent growth curve from experiments where larvae have been exposed to a relatively small number of constant temperature profiles. Since the temperature influences the developmental speed, a crucial step is the time alignment of the curves at different temperatures. We propose a model for time varying temperature profiles based on the local growth rate estimated from the experimental data. This allows us to estimate the most likely hatching time for a sample of larvae from the crime scene. Asymptotic properties are provided for the estimators of the growth curves and the hatching time. We explore via simulations the robustness of the method to errors in the estimated temperature profile. We also apply the methodology to data from two criminal cases from the United Kingdom.


Introduction
Forensic entomology is the study of insects (and other arthropods) in relation to criminal investigation, frequently involving insect evidence in cases of suspicious death. A great diversity of insects are attracted to decomposing human corpses, both to feed and to lay their eggs or larvae. Flies and beetles are the most common visitors, as both immature insects (larvae or maggots) and adults. In particular, blow flies are among the most important insects in criminal investigation, because they are usually the first to arrive (Greenberg, 1991). While the methods described in this work could be applied in general to any insects of forensic interest, we will consider here two species of blow fly, Calliphora vicina and Calliphora vomitoria (Diptera: Calliphoridae). Fly larvae are usually used to estimate a lower bound for the time since death, also known as post-mortem interval. Indeed, the age of the oldest insects on the body provides an assessment of when the mother insects gained access to the body and, as a consequence, a lower bound for the post-mortem interval in a criminal investigation. In this work we focus in particular on the information that is possible to obtain from the length of the larvae observed on the body or at the crime scene. The life cycle of flies is divided into different stages. In the case of an . Further details on the experimental procedure can be found in Donovan et al. (2006).
outdoor crime scene, adult female blow flies can arrive within just a few hours to lay eggs on the body (Reibe and Madea, 2010;Hofer et al., 2017). The eggs then hatch into first instar larvae (commonly called maggots) which start feeding and subsequently develop into second and then third instar larvae (between these stages, the cuticle is shed to allow for growth). When larvae have finished feeding, they usually disperse from the body and the cuticle shrinks and hardens, leading to the immobile puparia within which the pupal and pharate adult stages develop through metamorphosis . Finally, adult flies emerge from puparia, feed, mate and begin a new cycle. The rate of development of fly larvae is temperature dependent and when comparing the specimens observed at the crime scene with experimental data it is essential to adjust outcomes with respect to the different temperatures. Currently, the standard technique used to estimate the rate of development with temperature compensation is the Accumulated Degrees Hours (ADH) model. This consists in summing the hours of development multiplied by the temperature (in • C), thus providing a rough measure of the amount of thermal energy available to the larva. This can then be compared with existing experimental results about how many ADH are needed to reach each developmental stage in the insect life cycle and from this deduce the time of colonization at the scene (Amendt et al., 2007). Our goal is to provide more refined methods that consider the growth dynamic between hatching and pupariation, mostly ignored by previously existing techniques. We will rely on experimental development data on larval lengths for Calliphora vicina (see Figure 1) and Calliphora vomitoria (Richards et al., unpublished) to achieve this.
While experimental data are usually collected at constant temperatures, real life scenes are subjected to dynamic temperature profiles. In this work, we apply techniques from functional data analysis (FDA; Ramsay and Silverman, 2005;Ferraty and Vieu, 2006;Horváth and Kokoszka, 2012) to reconstruct the growth curve of the larvae in the presence of a varying temperature profile, using information from experiments run at constant temperatures (see again Figure 1). The objective is the assessment of the time since colonization, which is the interval between the earliest time when the observed blow fly specimens were laid as eggs and the time when they were collected at the scene. The current approach is often to use the ADH model which uses an average temperature approach. A more accurate assessment of the hatching time for a given set of larvae, via a case specific growth curve, would be a useful addition to the forensic science toolbox. This would be relevant for many criminal investigations, since evidence from forensic entomology is used in tens of high-profile cases every year in the United Kingdom alone, and many more worldwide (Tomberlin and Benbow, 2015). From a statistical perspective, FDA is the field that investigates the statistical properties of a collection of curves, surfaces or any more complex mathematical object. In forensics, FDA methods have already proved useful in analysing data generated from spectroscopic techniques (Dias et al., 2013;Burfield et al., 2015). More recently, FDA has also started to be of interest in forensic entomology, with one preliminary investigation using the technique for thermal wavelength analysis of maggot populations for time of death prediction (Warren et al., 2017). Here, we use the well known connections of FDA to growth curve analysis (Ramsay and Silverman, 2005). The main challenge in this paper is the need to combine different sources of information about the growth process to estimate the most likely growth curve that led to the larval lengths observed at the scene. Larval development data are collected in incubators which keep the temperature constant, while larvae at the crime scene are exposed to a varying temperature profile. The first part of the methodology requires an estimate of the actual case specific growth curve from the growth curves observed at constant temperatures in the lab. Then, the hatching time that leads to the best fit of the measured larval lengths at the crime scene is selected. The method has been developed with considerable attention paid to the pieces of information and data that are available in practice to forensic scientists. To support the validity of the proposed estimator, in addition to providing asymptotic results both for the estimator of the growth curve and of the hatching time, we assess the finite sample properties of the procedure via simulation studies. Finally, we discuss the application of the method to the data from two investigations in the United Kingdom. It is also worth mentioning that, while we developed our procedure in view of the application in forensic entomology which is described below, the methods and theory we discuss here for the estimation of growth curves can be relevant in other biological applications, such as microbiology (Zwietering et al., 1990) or quantitative genetics (Kingsolver et al., 2004).

Model and methods
The fact that larval development depends mainly on temperature is well known in forensic entomology (Donovan et al., 2006;Amendt et al., 2007). In addition, outdoor temperatures vary with time. However, laboratory experiments are bound to measure this relationship only for a relatively small set of temperatures, as shown in Figure 1. For each plot, the growth process is observed for a constant temperature profile (i.e. same temperature from hatching to pupariation). Our first goal is to use techniques from functional data analysis and nonparametric regression to estimate the expected growth curve corresponding to realistic temperature profiles. Once growth lengths are available for any temperature profile, one is able to estimate the hatching time and hence the interval between hatching and body discovery from data available at the scene.

Estimation of the growth process
Experimental larval development data. In a typical larval development data set, we have K experimental temperatures T 1 < T 2 < · · · < T K and, for each species of interest, one observes the larval lengths Y kjl measured at time t k j after hatching for j = 1, . . . , n k and for l = 1, . . . , N kj individual larvae at each time point which has been exposed to a constant experimental temperature T k , with k = 1, . . . , K. The observation times t k 1 , . . . , t k n k may differ across experimental temperatures, generally being set at wider intervals for lower temperatures and as a result of the experimental procedure a different number of individual larvae may be measured at each observation time (see for example Donovan et al. (2006) for more details). Figure 1 illustrates an example of this kind of experimental data for Calliphora vicina. In this example, we have K = 9 values of temperature set in the incubator, for each value of temperature T k we have a number of observation times n k which is between 16 (for T 6 = 18 • C) and 34 (for T 9 = 27 • C) and for most of observation times t k j we have about 15 measured larval lengths (fewer measurement are available for later times, as can be seen from Figure 1). In these data, we conventionally set the hatching time t h to be zero (i.e. t h = t 1 = 0) so that it is the time reference.
Estimating mean larval length curves for any experimental temperatures and their derivative. We can then assume that the observed lengths satisfy a nonparametric regression model Y kjl = L T k (t j ) + kjl , where kjl are independent, zero mean random variables and the mean larval length curve L T k depends on the experimental temperatures T k , k = 1, . . . , K. It is then possible to estimate L T k by means of a nonparametric smoothing estimator. In this work, we use a local linear regression estimator (see for instance Fan 1992, Fan and Gijbels 1992, Fan 1993, and Ruppert and Wand 1994. It is well known in the nonparametric statistical literature that local linear regression provides accurate estimates L T k for all growth curves observed at a constant temperature. However, as explained in the next paragraph, the introduced dynamic model for determining the hatching time requires the estimation of the derivative L T k of L T k and hence local linear regression is particularly suited to our setting. Figure 2 displays the estimated growth curves L T 1 , . . . , L T K at any time t in their respective time ranges for the Calliphora vicina data. Differential equation model. The question now is how to produce an estimate for the growth curve at a grid of time points t 1 , . . . , t M associated to a generic temperature profile T (t 1 ), . . . , T (t M ). We claim that, if we consider a small enough time interval, the temperature in that interval can be considered roughly constant and therefore the growth process would be bound to follow the local dynamics of the correspondent constant temperature growth curve at the corresponding stage of the growth process, i.e. at the point of the curve which reaches the current length in the growth process. This suggests the following differential equation model (DEM) for the local growth process: where L T (tm) is the growth profile at constant temperature T (t m ) and L(t) is the growth profile with varying temperature. The differential equation model (1) can be reformulated as L (t) = L T (tm) • L −1 T (tm) • L(t m ) where L stands for the derivative of L and the symbol • refers to the usual composition operator between real-valued functions. This dynamic model means that the (expected) local increment in length at time t m is the one that would occur in the growth profile Clearly the larvae at 4 • C developed with a different profile to those at all higher temperatures and none were observed to pupariate, most probably because 4 • C was near to the estimated lower development threshold for this species of 1.0 • C-1.5 • C (Donovan et al., 2006).
at constant temperature T (t m ) when the length is equal to L(t m ). This allows us to reconstruct the varying temperature growth profile iteratively, as a discretised solution to the above DEM: , and for m = 1, . . . , M − 1, with T (t 1 ), . . . T (t M ) being the varying temperature profile. This would solve the problem if we knew the expected growth curve L T and its derivative L T for any temperature T we can observe, but in practice we can have experimental data only for a relatively small set of temperatures.
Estimating growth profile and its derivative for any temperature. We need to estimate first the growth profile L T for a generic constant temperature T from a set of estimated growth curves L T 1 , . . . , L T K . The main difficulty here is that the temperature influences the speed of the growth process. Using the language of functional data analysis, these curves present both amplitude and phase variation (see Ramsay and Silverman, 2005;Marron et al., 2014), as can be appreciated from the example of the Calliphora vicina experimental data and their estimated growth curves at constant temperature (see Figure 2). At this stage, one postulates that the mean larval length curves L T 1 , . . . , L T K have corresponding profiles in a standardized time scale S T 1 , . . . , S T K . These can be obtained using warping transformations w T 1 , . . . , w T K acting over the time in such a way, for any k = 1, . . . , K and for all t in [0, t pup, where t pup,T k is the pupariation time (the last time in the experiment when the lengths of larval specimens were recorded) at constant temperature T k . The warping functions are assumed to be one-toone nondecreasing functions mapping [0, 1] into [0, 1] and restricted to the space of quadratic polynomials. They can be determined via a landmark registration procedure that aligns a set of structural points featuring the sample of curves (see, e.g., Ramsay and Silverman, 2005, Chapter 7). In our application, hatching time, time of maximum length and pupariation time are three structural points that correctly depict the growth profiles. We propose to match all growth profiles with respect to these three structural landmarks. The procedure achieving the alignment of a collection of profiles is well known in the functional data analysis literature under the terminology curves registration (see Kneip and Gasser, 1992;Ramsay and Li, 1998;Maldonado et al., 2002;Wang et al., 1997). It is of course possible to resort to more recent registration methods (see, e.g., James et al., 2007;Kneip and Ramsay, 2008;Srivastava et al., 2011) but the relatively simple shape of the growth curves allows landmark registration to provide easily and quickly two collections of simple curves: K estimated growth shapes S T 1 , . . . , S T K (i.e. the registered curves) and K estimated warping functions w T 1 , . . . , w T K (see Figure 3). From a practical point of view, the warping functions w T k 's are first computed (see Step 2 in Section 4 and the proof of Theorem 2 for their exact definitions) and then the growth shapes S T k 's are derived with the decomposition resulting from the registration procedure: for any u ∈ [0, 1], . The growth shape and warping function are reconstructed at any temperature T by functionon-scalar nonparametric regression (see Ferraty et al. 2011) as and for all u in [0, 1], where K S and K w are kernel functions and h S and h w suitably chosen bandwidths. The nonparametric kernel estimator is especially well adapted for the warping function setting. As a convex combination of the w T k 's, the estimator w T inherits the properties of the warping functions (i.e. w T is nondecreasing with w T (0) = 0, w T (1) = 1). The estimator of the growth process L T at any constant temperature T and for all t , where t pup,T is the predicted pupariation time depending on the temperature T . Given the sample of temperatures and corresponding pupariation time (T 1 , t pup,T 1 ), . . . , (T K , t pup,T K ), the pupariation time for a generic temperature T can be obtained with standard nonparametric regression as where K pup is a kernel function and h pup a data driven smoothing parameter. From now on, we are able to estimate the growth profile L T (t) for any temperature T and for all t in [0, t pup,T ]. According to the differential equation model (1), we also need to estimate the derivative L T of the growth curve L T . We know that for any t in [0, The previous steps provide w T , w T and t pup,T . The only quantity we need to recover L T is S T . For any u in [0, 1] and any k = 1, . . . , K, since the smoothing step provides L T 1 , . . . , L T K and the registration procedure gives the exact analytical writing of the warping functions w T 1 , . . . , w T K which enables the calculation of the inverse functions Given the sample of temperatures and corresponding growth shapes (T 1 , S T 1 ), . . . , (T K , S T K ), the computation of S T for any T can be derived from the function-onscalar nonparametric regression. Thus, one can deduce L T for any temperature T .
Conclusion and complement. The use of our DEM combined with the estimations of L T and L T for any temperature T enable the reconstruction of the growth curve for any temperature profile. The estimating procedure can be divided into four steps: (1) estimation of the mean larval length curves L T 1 , . . . , L T K and their derivatives L T 1 , . . . , L T K for each experimental temperature T 1 , . . . , T K , (2) decomposition of the mean larval length curves into growth shapes and warping functions, (3) estimation of the growth length L T and its derivative L T at any temperature T , (4) for a given temperature profile {T (t); t ∈ T }, computation of the corresponding length profile { L(t); t ∈ T } thanks to DEM (1). Asymptotic properties for all these estimators are described in Section 4 and stated with proofs in the Supplementary Material. However, a few technical issues need to be considered. First, as the growth curve is not monotone, there may be multiple times t at which L T (t k ) (t) = L(t k ). In particular, we need to distinguish between the feeding phase (initial monotone increase in length up to the maximum) and the post-feeding phase (the usually decreasing region after the maximum). To do this, as long as max l≤k L(t l ) < max t L T (t k ) (t) the process is labelled as increasing and the positive value of is used to obtain the length at t k+1 . Otherwise, the process is considered as post-feeding and the point corresponding to the negative derivative is used to update the length. Moreover, once the process is in the post-feeding phase, if the length reaches the minimum post-maximum value for the current temperature, i.e. if L(t k ) ≤ S T (t k ) (1), we assume that the process gets to pupariation and the larval length is not defined from that time on.

Estimation of the hatching time
We now want to use the estimated temperature-dependent growth profile to select the most likely hatching date given a set of length measurements taken at a time t * where the reference time is the local one. To do that, we compare the growth profiles that would be expected if the hatching time t h was at any time between the last time the victim has been seen alive t a and t * . Let L(t − t h ), t h ≤ t ≤ t * , be the growth curve for hatching time equal to t h and temperature profile {T (t * − t); t h ≤ t ≤ t * } . Let then Y * i , = 1, . . . , n obs be the measured larval lengths with i , i = 1, . . . , n obs independent random errors with zero mean and unknown variance σ 2 . Then, we can estimate t h as i.e. we choose the hatching time whose expected length at time t * best fits the observed values. However, we may also want to include some expert knowledge in the estimation procedure. First, the forensic entomologist collecting the sample may recognise whether the larvae reached the post-feeding phase, i.e. the region of the growth curve after the peak where larvae stop feeding in preparation for pupariation and subsequently decrease in length, or not. We can easily integrate this piece of information in the estimation procedure by restricting the admissible region for the minimisation problem (6) to the hatching times whose associated growth process at time t * has already had (or had not) reached the postfeeding region, i.e. the estimated derivative is negative at some time t ≤ t * . If we are willing to assume a parametric model for the error, we can also build an approximate confidence interval for the hatching time, for example by inverting the region of non rejection of the likelihood ratio test, i.e. the confidence region will be CR , where l denotes the log-likelihood and χ 2 1−α (1) the 1 − α quantile of the χ 2 distribution with 1 degree of freedom. However, this region is not necessarily convex and we may prefer a more conservative interval defined as CI(α) = [min CR(α), max CR (α)]. This is a connected interval that, by construction, asymptotically guarantees at least a (1 − α)% coverage of the hatching time. To be coherent with the estimation procedure above, we also need to account for the information about the developmental stage of the observed larvae. Since we optimised the criterion only in the admissible region of hatching time which guarantees the correct developmental stage (postfeeding or not) at the time of sample collection, we need to do the same for the parameter space where the log-likelihood is defined. This is then restricted to the same set of hatching times used in the estimation. On the other hand, one may want to include prior information about the hatching time, coming for example from the investigative activity (in addition to the interval of admissible hatching times). Let us assume we can translate this information into a prior distribution on the parameter t h , so that t h ∼ π where π is a known distribution. We can then use Bayes theorem to update the information about t h given the observed larval lengths and derive a posterior distribution for t h . For example, assuming a normal distribution for the errors we have In practice we need to substitute σ 2 with a plug-in estimate, for example the sample variance of the observed lengths and L(t * − t h ) with the estimated quantity L(t * − t h ). Note that in the estimation procedure in this section, we are ignoring the uncertainty in the estimation of the time-dependent growth process. This is indeed present and, while the uncertainty from the experimental data is usually negligible, errors on the temperature at the crime scene can affect the results. We are going to explore this issue through simulation studies in Section 3.1. Note that it is straightforward to generalise these ideas to the case where more than one species of larvae (for which developmental data are available) are observed. Let Y * ij be the observed length of the i-th sample from the j-th species, j = 1, . . . , J, with ij , i = 1, . . . , n j independent random errors with zero mean and unknown variance σ 2 j . Then, where σ 2 j is the sample variance of the lengths of the j-th species. To support the validity of the proposed estimator in practice, we investigate its properties in three different ways. Simulation studies and applications of the method to a couple of police case studies can be found in Section 3, while Section 4 is devoted to the exploration of the theoretical properties of the estimator in the asymptotic regime. In particular, it is shown in Theorem 4 that the rate of convergence of the estimated hatching time t h especially depends on the number of larval lengths observed at the crime scene, the grid size where the temperature profile at the crime scene is sampled and the number of experimental temperatures. These quantities have to be as large as possible to ensure an accurate estimation of t h .

Empirical demonstrations
In this section, we first consider some simulation studies to assess the robustness of the proposed method to errors in the temperature profile measured at the crime scene. Then, we present the application of the method to the data coming from two investigations. For reasons of privacy, the cases have been anonymised and we conventionally set the time at which the measurements of larval lengths have been taken to be t * = 0, with negative times indicating the hours before this moment. In both these forensic cases, we do not have an external corroboration (such as a defendant confession) of the time the body has been abandoned, therefore we compare the results of our procedure with those provided by the ADH method currently used in criminal cases. Throughout this section, we use quartic kernels in the estimation of the growth shape, the warping function and the pupation time, with the same bandwidth h = 4 • C. This is chosen by visual comparison between the estimated growth curves and the observed experimental curves, as it can be seen in Figure 8 and 10 for the growth profiles of the real criminal case studies.

Simulation studies
We illustrate here the performance of the proposed methods. We first simulate from model (5) using two different temperature profiles. The first, scenario (a), is a constant temperature profile set to 10 • C (i.e for all t, T (t) ≡ 10). The second, scenario (b), corresponds to a more realistic varying temperature profile which is a subset of the temperature measurements from the weather station from the second case study in Section 3.3. In both cases, an hourly time grid t 1 = −200, . . . , t 201 = 0 is used to evaluate the growth curves and optimise the criterion (6) where true hatching time t h is set to −100. We also assess how the results are sensitive to variations in temperatures. To this end, noisy temperatures are build with an independent Gaussian error: . We consider four different values for σ T = 0.1, 0.25, 0.75, 1. For each value of σ T , we simulate 1000 samples of 20 observed lengths and we estimate 1000 hatching times. Figure  4 give an idea on the shape of temperature profiles under scenario (a) and (b), respectively, and of the corresponding estimated temperature-dependent growth curve. Figures 5 and 6 show the distributions of the estimated hatching times for the different values of σ T . The dispersion of the estimates is smaller for scenario (a) which is in accordance with the basic structure of the constant temperature profile. The variability of the estimates also increases with the variance of the error added to the temperature (i.e σ T ). We then compare our proposed methods with existing techniques in these simulated scenarios. While the most commonly used approach based on the accumulated degrees hours (ADH) model is usually carried out by leaving the larvae in an incubator until pupation, we could use as benchmark a relationship between ADH and length estimated from the experimental data. We use experimental data on Calliphora vomitoria growth in Richards et al. (unpublished) to estimate this relationship nonparametrically via smoothing splines using the gam function from the R package mgcv. Figure 7(b) shows the raw data and the estimated curve. We then use this estimated curve to predict the hatching time associated with the ADH predicted in correspondence of the average observed larval length. We set a larval development threshold (i.e. the temperature below which there is no growth in the larvae) to 4.3 degrees Celsius, as observed experimentally for the larval development phase up to pupariation. We consider here only the more realistic scenario (b) of varying temperature profile and the results are shown in Figure 7(a). The proposed method based on the reconstruction of the growth curve appears to be more accurate than the one based on the ADH curve. Although the ADH estimation can be better tuned by changing the lower developmental threshold (and potentially having it non constant along the larval development), this is often challenging in practice and it would require to develop a model for ADH similar to the one we propose here for the growth curve.

First case study
We consider a first case study where n obs = 70 Calliphora vicina post-feeding larvae were collected from the body. In this investigation, there was not a unique crime scene since the body was moved between death and discovery. For this reason, the temperature profile to which the body was subjected is provided here by forensic experts based on the information about the body location coming from the investigation. Figure 8 illustrates the temperature time series for the 371 hours before the time the larvae were killed, prior to subsequent measurement, together with the constant temperature growth profiles corresponding to each observed temperature in the series. The interval of 371 hours was the largest one that was considered possible by forensic scientists at the scene. For this case, the application of the ADH method suggested as plausible interval for eggs' hatching the one between 275.5 and 208.5 hours before the measurements were taken (Donovan et al., 2006). Figure 9 shows the growth curves for the hatching time estimated with (6) and the profile of the objective function in the minimisation problem. The estimated hatching time is -266 hours (before the larval measurement), which is within the range obtained from the ADH method. However, the flat plateau in the criterion suggests little stability for the estimate. Indeed, if we assume a Gaussian distribution for the measurement errors, the 95% approximated confidence interval procedure gives us an interval of [−304, −193], which includes the range suggested by the ADH method. Note that here the computation of the criterion is restricted to the admissible region of hatching times for which the expected growth curve would have reached the post-feeding stage by the time the larvae were collected.
In conclusion, the estimate for the hatching time provided by the proposed method is in agreement with the one from the ADH method in this case study, but the estimated uncertainty is larger than the bounds provided by accepted procedures based on ADH, for which no rigorous assessment of the uncertainty is available. Also, the fact that the point estimation of hatching time based on the proposed method gives us a larger time interval than the one obtained from  Figure 8: Left: Temperature profile to which the body was subjected in the hours before larval lengths were measured. Note that between the body discovery and the collection of the larvae the body was stored in a fridge for a few hours. Right: Estimated constant-temperature growth curves for each temperature in the observed interval (black lines), the temperature corresponding to observed experimental temperature are highlighted. Note the constant (zero) growth corresponding to the temperature in the fridge.
the ADH method is consistent with the simulation results in Figure 7, where the ADH method underestimated this time interval.

Second case study
This second case study also includes the estimation of the temperature profile. A logger measures the temperature at the crime scene after the body is discovered and this is compared to the data from the closest weather station, which are then used to estimate the past temperature at the crime scene. In particular, we are going to use a local polynomial kernel regression to estimate the crime scene temperature. The temperature at the crime scene, the weather station temperature measurements and the estimated crime scene temperature can be seen in Figure  10, together with the constant temperature growth curves for all the estimated temperature values. At the crime scene, the lengths of 9 Calliphora vomitoria larvae still in the feeding phase were measured (among specimens of other species, for which experimental developmental data are not currently available). We are using here the development data on Calliphora vomitoria collected at the Natural History Museum, London (Richards et al., unpublished). The assessment from forensic scientists, based on the application of the ADH method as well as qualitative considerations about the other species present at the scene, is that the body had been infested by bowfly eggs between 270 and 240 hours before the body was discovered. Figure 11 shows the estimated growth curve and the profile of the objective function. The estimate for the hatching time is −251 hours before the body discovery and the 95% approximate confidence interval for a Gaussian error model is [−252, −247]. However, note that −252 is on the boundary of the admissible region and this makes the confidence interval based on the likelihood ratio statistics invalid. Moreover, we can see that the final expected length is too low with respect to the data and indeed the minimum for the criterion is reached at the boundary of the admissible region. This may suggest that the actual temperature at the crime scene was higher than predicted or that larvae reared in experimental cultures do not reach the potential maximum size reached in natural conditions. This is, of course, valuable information for the forensic scientists, unavailable from the ADH model. The original forensic reports proposed a time immediately adjacent to the time of the victim disappearance. This was because information from the investigation about the last time the victim has been seen alive was also taken into account. We can do the same with our model by using a Bayesian approach when the prior is taken to be a uniform between the earliest possible hatching time after the last sighting of the victim alive and the time of discovery of the body. The results can be seen in in Figure 12. Our method suggests strongly the earliest possible time for the hatching, in agreement with experts' judgement. It should be noted again that the growth curve does not appear to be able to reach the observed lengths. In conclusion, the experimental population did not develop to the average size of the larvae observed at the crime scene. This may be either because the development of larvae at the crime scene were affected by something that was not considered in the laboratory experiments (for example, larval-generated heat due to high concentration of specimens) or because the temperature at the crime scene was higher then expected during a portion of the development period. However, one useful advantage of the method over the ADH approach is that the estimated growth curve provides a diagnostic tool to spot any issue with the assumptions on the growth process.

Asymptotic properties
The main goal of this section is to state the asymptotic properties of the estimated hatching time, giving theoretical justification for our estimators, in addition to the empirical assessment given in the previous section. To this end, the estimating procedure is decomposed into four steps. In the first step the constant temperature growth profile L T k (·) and its derivative L T k (·) are estimated from observed growth data set at a given ambient constant temperature T k belonging to the fixed design T 1 , . . . , T K . A decomposition of the estimated constant temperature growth profile L T k (·) into growth shape S T k (·) and warping function w T k (·) is achieved in step 2, at a given ambient temperature T k taken in the same fixed design.
Step 3 focuses on the estimation of the constant temperature growth-profiles L T (·) and its derivative L T (·) at any ambient temperature T (i.e. T may be outside the experimental design T 1 , . . . , T K ). In the last step the asymptotic behaviour of the estimated hatching date t h is provided. Assumptions are formally stated in the appendix and proofs can be found in the Supplementary Material. In order to make the theoretical developments more readable and accessible, the framework is voluntarily restricted to fixed designs, non random temperatures and only one species of fly. However, the setting can, of course, be extended to random designs, random temperatures and several species of fly, subject to some straightforward adjustments. Because all regularity assumptions involving L T , S T , and w T are uniform over a compact set of temperatures, all asymptotic results are also uniform with respect to the temperatures.
Step 1: Estimating constant temperature growth profiles. For a given constant temperature T k ∈ {T 1 , T 2 , . . . , T K } one observes repeated larvae growth lengths {Y kjl ; j = 1, . . . , n k , l = 1, . . . , N kj } at a grid of n k times after hatching 0 = t k 1 < t k 2 < · · · < t k n k := t pup,T k . From this sample, one derives the growth profile L T k (t) at constant temperature T k and time after hatching t by means of the nonparametric regression model Y kjl = L T k (t k j ) + ε kjl and the local linear regression estimator ( a, b) = arg min with where, for any k, Y kj : , and the (n k × n k ) diagonal matrix . In practice, the function ω k is set in such a way ω k (t k j ) = N kj / j N kj . The expression of the estimator of L T k (t), the first derivative of L T k (t), is given by: As it is shown in the next section, to derive the consistency for the whole estimating procedure, the uniform rate of convergence for L T k and its derivative L T k are needed.
Theorem 1 Set r n = Under (H1), (H4)-(H8), for any T k ∈ {T 1 , . . . , T K } it holds: and Although the literature on local linear regression is dense, the proof of (10) and (11) needs some adjustments that are presented in the Supplementary Material.
Step 2: Decomposing growth profiles into growth shape and warping functions.
For the growth profile L T k (·) at a given constant temperature T k , one postulates the existence and unicity of a growth shape S T k mapping [0, 1] into R and a strictly increasing warping function w T k mapping [0, 1] into [0, 1] such that, for any time t after hatching, L T k (t) = S T k • w T k (t/t pup,T ), where S T k and w T k are derived from L T k (·) by aligning the curves at 0, t k max /t pup,T and 1. For a given real α ∈ (0, 1) and for any k varying from 1 to K, the warping function w T k is assumed to be a quadratic polynomial such that: In other words, the warping function w T k is a strictly increasing function interpolating the three points with coordinates (0; 0), (t k max /t pup,T ; α), and (1; 1). Imposing w T k in the space of quadratic polynomials ensures its unicity. This procedure is commonly called registration (see Ramsay and Silverman, 2005). So, our aim is to state asymptotic properties of the K estimated pairs ( S T 1 , w T 1 ), . . . , ( S T K , w T K ) and their corresponding derivatives ( S T 1 , w T 1 ), . . . , ( S T K , w T K ). For any k ∈ {1, . . . , K}, let t k max be the time of maximum length derived from the estimated growth profile L T k ; the estimated warping function w T k is defined as the strictly increasing quadratic polynomial satisfying w T k (0) = 0, w T k ( t k max /t pup,T ) = α, w T k (1) = 1, and, for any u in [0, 1], the corresponding estimated growth shape Theorem 2 Under (H1)-(H8), for any T k ∈ {T 1 , . . . , T K }, w T k (resp. S T k , w T k , S T k ) converges to w T k (resp. S T k , w T k , S T k ) with the same rate of convergence r n = O(h 2 L ) + O P (nh 6 L ) −1/3 .
Step 3: Estimating constant temperature growth profiles and corresponding first derivatives at any ambient temperature. In this step, the challenge is to estimate the constant temperature growth profile L T at any ambient temperature T (i.e. not only in the discrete grid T 1 , . . . , T K ). Based on the K-samples ( S T 1 , T 1 ), . . . , ( S T K , T K ) and ( w T 1 , T 1 ), . . . , ( w T K , T K ) one is able to derive an estimator S T (resp. w T ) of the growth shape S T (resp. w T ) for any ambient temperature T (see (2) and (3)). So, for any temperature T , it is easy to deduce an estimator of the growth profile L T by setting L T (t) = S T • w T (t/t pup,T ) for any t in [0, t pup,T ]. Similarly, based on the K-sample ( S T 1 , T 1 ), . . . , ( S T K , T K ) and ( w T 1 , T 1 ), . . . , ( w T K , T K ), one can derive the estimation of S T and w T at any temperature T by setting where K S (resp. K w ) is also a kernel function and h S (resp. h w ) the non-negative smoothing parameter. Because for any temperature T , the derivative of the growth profile L T (t) for any t in [0, t pup,T ] is estimated with L T (t) = S T w T (t/ t pup,T ) w T (t/ t pup,T ) t −1 pup,T .
Before tackling the estimation of the unknown hatching date t h , one focuses on the varying temperature growth length L(t) at the time t after hatching which is assumed to depend on the temperature process {T (v), v ∈ [0, t)}. One postulates that the growth length L(t) given the set of temperature variations {T (v), v ∈ [0, t)} satisfies the dynamic growth model: where L(0) = L T (0) (0) is the length at the hatching time. This formulation can be expressed in terms of differential equation: In practice, one uses a fine grid of time 0 = t 0 < t 1 < t 2 < · · · < t M −1 < t M < t < t M +1 such that sup |t +1 − t | = O(M −1 ) and the previous dynamic growth model can be approximated by its discretized version: with L(t 1 ) = L T (t 1 ) (t 1 ). So, the estimator L of the varying temperature length profile L is defined as follows: Once the varying temperature growth profile is estimated, the ultimate task is to compute the hatching time and the hatching interval H. To this end, one observes on the crime scene at a given date t * a sample of i.i.d. larval lengths Y * 1 , . . . , Y * n obs corresponding to the growth length L(t * − t h ) reached at time t * − t h (:= H) after hatching and computed from the past outdoor temperature time series {T (t * − t); t ∈ [t h , t * ]}, where t h stands for the unknown hatching date. Moreover, for i = 1, . . . , n obs , one assumes Y * i = L(t * − t h ) + i where the i 's are zero mean random errors. Then, the estimated date of hatching t h is the minimizer of , which gives the estimation of the hatching interval: Theorem 4 Under (H1)-(H8) it holds: As a by-product, one gets under same assumptions the consistency of the hatching interval:

Conclusions
We have described a functional data approach to incorporate the information from the constant temperatures used to gather experimental data into the estimation of the crime scene varying temperature growth profiles. This can be used to estimate the most likely hatching time for the observed larvae. The proposed method has some advantages over the existing ADH approach. First of all, it can be applied directly to the lengths of the larvae measured at the crime scene, without the need to wait for the larvae to develop up to the next stage of the life cycle in incubators, as can be required by ADH methods based on duration of discrete life stages. Second, the proposed method allows forensic scientists to also consider the estimated growth curve as a model diagnostic tool to highlight any problematic situations, as seen in the second case study. Moreover, when larvae are observed in the middle of the development process, the method provides a more accurate estimate of the hatching time together with an estimate of the uncertainty. Note that this was not the case in the two case studies we considered, where larvae had already reached the maximum size or were in the post-feeding phase. We have demonstrated both theoretically and empirically that the functional data approach can provide good estimates for the interval from hatching time to discovery of the body. These are, of course, subject to a number of uncertainties, such as those relating to temperature estimation as well as the inherent biological variability in both the experimental and crime scene entomological data. We have shown that this uncertainty can be captured both in a frequentist or Bayesian setting, allowing the model to be used in either framework, depending on the availability of additional prior knowledge and the legal requirements of the court. It needs to be stressed that the methods by no means give a definitive conclusion about the post mortem interval on their own, since, on the one hand, many factors can delay the access of flies to the body (see, e.g., Bhadra et al., 2014) and, on the other hand, only expert judgement about the surroundings can guarantee that the observed larval specimens are the oldest to have colonised the body. However, they do yield more complete estimates of the growth of the larvae compared with the simpler accumulated degree hour models, and therefore are of use in fieldwork situations. For the future, we also need to expand the model to include the egg and intra-puparial stages of immature blow fly development, so that the proposed methods extend up to the time of colonisation. Further work would also consider a more comprehensive handling of the uncertainty in the Bayesian framework to account for the uncertainty in the estimation of the average growth process conditional on the temperature profile, as well as potentially the uncertainty on the temperature profile itself. Finally, the selection of the optimal bandwidth in the non-parametric regression is also still an open problem, this potentially being different for shape and warping functions and pupariation time.
About the sample sizes, grid size and the bandwidths used in the estimating procedure, one requests: (H7) set n := inf k n k ; n tends to infinity, h L tends to zero with n, and n h 6 L tends to infinity with n, (H8) the number K of experimental temperatures tends to infinity with n; h pup , h S , h S , h w , h w tends to zero with K and K h pup tends to infinity with K; the size p of the grid where the generic temperature profile is sampled tends to infinity; the sample size n obs of larval lengths observed at the crime scene tends to infinity.