Sensitivity of temperature-based time since death estimation on measurement location

Rectal temperature measurement (RTM) from crime scenes is an important parameter for temperature-based time of death estimation (TDE). Various influential variables exist in TDE methods like the uncertainty in thermal and environmental parameters. Although RTM depends in particular on the location of measurement position, this relationship has never been investigated separately. The presented study fills this gap using Finite Element (FE) simulations of body cooling. A manually meshed coarse human FE model and an FE geometry model developed from the CT scan of a male corpse are used for TDE sensitivity analysis. The coarse model is considered with and without a support structure of moist soil. As there is no clear definition of ideal rectal temperature measurement location for TDE, possible variations in RTM location (RTML) are considered based on anatomy and forensic practice. The maximum variation of TDE caused by RTML changes is investigated via FE simulation. Moreover, the influence of ambient temperature, of FE model change and of the models positioning on a wet soil underground are also discussed. As a general outcome, we notice that maximum TDE deviations of up to ca. 2–3 h due to RTML deviations have to be expected. The direction of maximum influence of RTML change on TDE generally was on the line caudal to cranial. Supplementary information The online version contains supplementary material available at 10.1007/s00414-023-03040-y.


Introduction
Temperature based death time estimation (TDE) is crucial in homicide investigations.TDE is estimated from core temperature measurement data using the model curve T(t) comprising a priori knowledge of the postmortem rectal temperature decline.The phenomenological approach and the physics-based approach are two different techniques used to generate the model curve T(t).The method of Marshall and Hoare with the Henßge parameters from the Nomogram method (MHH) is a prominent phenomenological approach using a double exponential model with fitted parameters for modeling rectal cooling [1].Physics-based approaches use the heat transfer equation considering heat exchange mechanisms, cooling conditions and thermal material properties.Finite element (FE) based TDE method (FEM) is a physicsbased approach with reasonable computational effort (see e.g.[2][3][4]).For an extensive overview over temperature based TDE see e.g.[5].
In any TDE method, error quantification is a desideratum as TDE results can lead to acquittal or conviction of suspects.Three different errors exist in FEM such as errors due to space and time discretization, input data errors, and model errors.They lead to the deviation of the TDE value from the actual value.Though several approaches for error estimation exist for TDE e.g. for MHH [6] and for FEM [7] none of them takes into account the uncertainty introduced in TDE due to input data errors caused by variations in the rectal temperature measurement location.The present article tries to close this gap using an FE approach.This does not, however, limit the results to FE TDE methods.
Rectal temperature measurement (RTM) can be performed using approved devices and without the necessity to injure the body.However, there is a certain amount of diversity in the characterization of the measurement locus in forensic literature.In early well-known TDE studies [8,9], it is reported a 3 -4 inches [7 -10 cm] insertion depth for RTM.Marshall and Hoare [10] mention the relationship between the liver and rectum temperature.Henßge [11,12] states that the insertion depth for RTM should be at least 8 cm from the musculus sphincter ani.In [11] Henßge estimates the precision of rectal temperature measurement to be about " + -2 °C".Further, Henßge [13] emphasized the usage of linear, rigid and not very flexible thermometers to achieve an insertion depth without applying force.He also recommends a temperature measurement at the mesenteric root in the lower left abdomen additional to the temperature measurement in the deep abdominal space [14].Moreover, Henßge advises the thermometer to be inserted in the rectum as deep as possible without applying any force [15].Indeed, a lot of discrepancy in rectal temperature measurement location exists in literature, so it is necessary to study its effects on TDE.Hence, this motivates a sensitivity analysis of TDE with respect to measurement locus.
The study was executed with participation of the Institute of Legal Medicine (IRM) at the University Hospital Jena of the Friedrich-Schiller-University Jena and of the Zuse Institute Berlin (ZIB), an interdisciplinary research institute for applied mathematics and data-intensive high-performance computing.
Additional computation results are presented in the Supplementary Information (SI) associated with the electronic edition of the article available via the website of the journal.

Finite element model
Background information on numerically solving partial differential equations and on FE simulation can be found e.g. in [16].
For one exemplary case of human corpse with a given CT-scan, two different methods were used to generate FE models.In one of the methods, the FE Model (CTM) was developed from the segmentation of a CT scanned human body and in the other, the human body was manually approximated by hexahedral elements (see e.g.[2,17]).Two variants of this manually generated model were constructed: One model (CM) is floating freely in air and a copy (CMS) of CM is laying firmly on a wet soil substrate.Detailed simulations of physical heat transfer processes generated the rectal temperature curves T(t).The corpse cooling was computed in all of the FE models based on the well known heat transfer Eq. ( 1), a partial differential equation where c is the specific heat capacity, ρ the mass density, and κ the heat conductivity of the tissue: Heat transfer from the body to the environment across the skin due to convection and surface-to-ambient radiation was captured by a Robin boundary condition with effective heat transfer coefficient γ as given in Eqs.(2), (3) (see e.g.[7]): Here, T A is the ambient temperature, h is the heat transfer coefficient, ε is the emissivity and σ is the Stefan-Boltzmann constant.In the CTM developed from the CT scan of a human body, the initial temperature field T 0 defined at time t = 0 satisfies Pennes' Bio-Heat-Transfer-Equation (BHTE) [18] where ρ b c b is the heat capacity of blood, w is tissue perfusion, and T core is the body core temperature.
In the manually developed FE model CM, the initial temperature field T 0 is defined with a gradient between core and outside elements as in [2] referring to physiology literature.
For reliable comparison between the different models, the equivalent effective heat transfer coefficient was applied to the CTM corresponding to the convection and radiation terms applied on CMS and CM.The model curves T(t) were sampled at the designated nodes C and SP R k (k = 1,…6) corresponding to the anatomical RTM locations (RTML).The point C was placed at the intended measurement location in the respective FE model, whereas the six points SP R k (k = 1,…,6) lay on the vertices of an octahedron with radius R. The SP R k were defined to represent possible locations of a misplaced temperature sensor.If no confusion can arise about the value of R the SP R k are abbreviated SPi.
Our results consist of maximum deviations D MAX M comparing two cooling curves computed at two locations of RTM's in a specified FE model M = CTM, CM, CMS.The distance D MAX M is defined between the cooling curves T C (t) and T SPi (t) at the center point C of a measurement point octahedron (see Fig. 4) of radius R and at its vertices SP R k (k = 1, …, 6).We also computed differences D MAX M1,M2 : = D MAX M1 -D MAX M2 of those maximum deviations for two pairs (CM-CMS and CM-CTM) of different FE models (M1, M2) for all C -SP R i pairs.All computations were performed for different constant ambient temperatures T A and for different octahedra radii R. We also performed all of the computations for each of four Q-ranges, where the first three are Henßges Q-ranges (see e.g.[15]).

CT meshed FE model CTM
All of the modeling work done on the CT meshed FE model CTM was performed at the Zuse Institute Berlin (ZIB).The CT scan of a male corpse of length L = 1.74 m and weight M = 62 kg was segmented in the software AMIRA and the segmented data was converted into an FE model, comprising of 961,234 tetrahedral elements, via the FE-program Kaskade (see e.g.[19]), an in-house code of ZIB.As depicted in Fig. 1, the feet of the corpse were not included in the model and the weight, size and volume of the whole model were scaled according to the scanned corpse dimensions.The segmentation of different body parts was carried out using differences in density and location of different tissue types like bone, fat, bladder, kidney, abdomen, liver, heart, lungs and muscle as seen in Table 1 in the Appendix..In CTM, skin tissue was not segmented separately from the subcutaneous fat tissue since skin tissue was assumed to have a negligible effect on the TDE.Hence, the skin material properties were not considered in CTM.

Coarse meshed FE model CM
Our coarse meshed FE model CM was already described in several publications [2,3,20,21].Its mesh was manually constructed using the FE preprocessor Mentat.It consists of approx.12,000 nodes in approx.8300 hexahedral elements with trilinear shape functions.The model contains numerous compartments (see Fig. 2) standing for distinct anatomical structures of the body.
Our standard FE model of length L = 1.64 m and mass M = 64 kg was scaled geometrically to a length of L' = 1.74 m and weight of M' = 62 kg, corresponding to the male reference body of our study.The scaling factors are defined as in the Eqs.( 5) and (6) where k 1 is linear 1-dim scaling along the body length axis and k 2 is linear 2-dim dilation in the transverse plane [2].
Generation and scaling of FE mesh and post-processing were performed in MSC's pre-and post-processing tool Mentat.
Further, to investigate the TDE sensitivity on more complex boundary conditions favoring thermal energy transfer from the body core in dorsal direction, a supporting structure or floor is modeled with an hexagonal FE mesh as shown in Fig. 2a.The FE model lies on its back on the support structure defined, presumably causing large differences between C-temperature and temperature at SP4 or SP3, i.e. along y direction.The support structure is defined with the thermal properties of wet soil such as thermal conductivity c = 2 W/ m 2 K, specific heat κ = 2200 J/kgK, density ρ = 1900 kg/m 3 [22] and emissivity ε = 0.95.The coarse model laying on its back on the support structure is abbreviated CMS, whereas the coarse model without support structure is named CM ( 5) (see Fig. 2c).CMS and CM were scaled to the same height and weight as CTM.This will eliminate the influence of weight on the body cooling.

Definition of measurement positions C and SP R k
TDE estimations are usually based on a single RTM.Hence, the measurement position is subject to variations due to different anatomies, thermometer angles and insertion depth.An approximate ideal measurement position can be specified from anatomy including forensic knowledge on measurement locations used in practice.The central position C for our measurements was determined in the pars ampullaris of the rectum near the incisura transversalis, which is the passage from pars ampullaris to the pars sacralis of the rectum (Fig. 3).
Additional measurement positions SP R k (k = 1, …,6) for our sensitivity studies were chosen in spatial x, y and z direction at distances (octahedral radii) R = 0.5 cm, 1 cm, 2 cm from the central measurement position C.For each fixed distance R six additional measurement locations SP R k were established.The six additional measurement positions SP R k form an octahedron with the central point C in its middle.The schematic sketch representing the additional positions in the model was illustrated in Fig. 4.
The anatomical position of additional measurement points in CTM are described as follows: SP1 lies in the center of abdomen surrounded by abdomen tissues, SP3 lies nearby bone, SP4 lies in the center of abdomen, SP6 lies very near to abdomen. Figure 1b shows the anatomical positions of measurement points in CTM.Similarly, the anatomical position of additional measurement points in CMS and CM are depicted in Fig. 2c, which helps to identify the position of additional measurement points.SP1 lies lateral from C in positive X direction.Due to CM's left-right symmetry the temperature value at SP2 is equal to the temperature at SP1.The points SP2, SP3, SP5 lie in the gastrointestinal tract.SP6 and SP4 lie on the border between gastrointestinal tract and pelvis.

Simulation
The CT-generated model CTM was imported from the FE research code Kaskade 7 [19] into the commercial MSC-Marc Mentat FE system and the material parameters are defined as given in Table 1 of [2] (see Table 1 in the Appendix).The same material properties were used also for CMS and CM models.The initial temperature field T 0 (r) at the time of death on the body was calculated via the Pennes' Bio-Heat-Transfer-Equation BHTE on the original CTM model [18].The field T 0 (r) was approximated by 10 discrete node sets corresponding to 10 discrete temperatures in the converted version of CTM.In the CMS and CM models, the initial temperature was defined .Although the method of defining boundary conditions differs between CTM and CM models, it was found that the linearization of the radiation term applied in CTM has negligible effect on the cooling curve [7].A

Quantification of TDE deviations
TDE deviation from the central position C to any position SP R k is quantified by the maximum distance D MAX between the respective cooling curves, which is defined in this section.Let the time interval [a, b] be given during body cooling which starts at 0 h and let a = 1 h be the left point and b = 45 h be the right point of our interval.Let further T 1 (t) and T 2 (t) be two cooling curves defined on the interval [0, b].Evaluation of the cooling curve distances started at a > 0 h for reasons of numerical stability.The ambient temperature T A is assumed constant in space and time outside the body.The reference curve T R (t) for two cooling curves T 1 (t) and T 2 (t) is constructed as the pointwise mean of T 1 (t) and T 2 (t): Let T 0 : = T R (0) be the initial temperature at the center point C. Now the reference curve T R (t) can be normalized to the function Q(t) taking values in the real number's interval [0, 1]: This makes it possible to define the Q-intervals Q 1 to Q 4 (where Q 1 , Q 2 , Q 3 correspond to Henßge's [15] normed temperature ranges for tolerance radii) and their left and right boundaries q 0 , …, q 4 by: Assuming T R (t) to be strictly monotonically decreasing, we can uniquely map the Q-interval boundaries q n to the respective temperature values T R(4-n) .
Since the two cooling curves T 1 and T 2 are monotonically decreasing it is possible to define their inverse curves t 1 , t 2 , called the (absolute) time since death (TSD) estimation result curves w.r.t. the cooling curves T 1 , T 2 on the domain [T MIN ,T MAX ] which is the range of the reference curve T R on the temperature axis by: The time functions t 1 , t 2 are frequently concatenated with the reference curve T R , thus giving t 1 (T R (t)), t 2 (T R (t)) in the following definitions.For short those concatenations will be abbreviated t 1 (t), t 2 (t).If the point t in time is clear, we will even write t 1 , t 2 .We will now quantify the distance between the cooling curves T 1 and T 2 , which are actually constructed as distances between the inverse curves t 1 and t 2 .The reason for this construction lies in our primal interest in time differences because our research's number one target is time since death.All of the temperature curve distance measures are constructed using the (time oriented!)maximum distance D MAX (T 1 ,T 2 ) between two cooling curves T 1 and T 2 , which is shown here for continuous functions first.The maximum domain of the reference curve T R on the time axis is the interval [t MIN ,t MAX ] which is taken therefore as the domain for maximum-finding in our definition: We will now consider the equivalent of the definition (14) in terms of real measured temperature curves T i (t).Each curve T i (t) is represented as a finite series of real measurement values (T i 1 , …, T i N ) on a finite (regular) grid of time values (t 1 , …, t N ) which is the same for both of the curves.The first point t 1 of the time grid is identical to the starting time of cooling computation, while the last point t N marks the end of the cooling computation interval.To provide a joint domain of definition for both our inverted temperature curves T R −1 we define: Now we have constructed the range of our reference curve T R (t).The limits of its domain of definition [t MIN ,t MAX ] were computed: Let K be a natural number and let (t 1 , …, t K ) be a regular time grid on the interval [t MIN ,t MAX ] with t 1 : = t MIN and t K : = t MAX .This provides the equidistant sampling points for the final definition of the functions t i (T R (t)).Let further be (T 1 , …, T K ) the corresponding temperature grid with T 1 ≤ … ≤ T K and T k : = T R (t k ) for all k = 1, …, K and T K : = T MAX .The time grid's width Δt is: Let further be t qn the inverse image of T Rn under the reference function T R for all n = 0, …, 4: The five points t qn on the time scale constitute the endpoints of four time intervals t Q1 , …, t Q4 : We will redefine our measure (14) quantifying the distance between the cooling curves T 1 and T 2 in terms of the real samples (t 1 , …, t K ) and (T i 1 , …, T i K ).The distance measure is the distance D MAX (T The (absolute global) maximum distance D MAX (T 1 , T 2 ) is defined as:

Results
In a first step we evaluated the results for the direction in space (X,Y,Z as in Fig. 4), where the largest TDE-deviation value D MAX caused by RTML -variation was seen.Interestingly, this maximum direction depended neither on the ambient temperature T A = 5 °C, 15 °C, 25 °C nor on the octahedron radius R = 0.5 cm, 1 cm, 2 cm.The FE-model chosen makes the only difference.Therefore we get the following list of directions for maximum D MAX caused by RTML-variation:  Results concerning relative TDE difference distances and distances based on the L 2 -norm instead of the MAX norm are presented in the electronic SI.

Discussion
Multiple analysis was carried out to determine the TDE error caused by measurement locus variations for the different FE-model types such as CTM, CM, CMS using maximum distance D MAX .
Our simulation results show maximum TDE deviations D MAX with respect to variations in RTML up to 2 cm in the order of magnitude of 2 to 3 h.Comparing to other influences like variations in ambient temperature T A (see [6,7]), initial body core temperature T 0 (see [20]), etc. the TDE deviations of 2 to 3 h are not negligible and can cause considerable errors in TDE.Our RTML-caused TDE deviation results D MAX lie within the 95% confidence interval of Henßge [15].For the first Q interval, Henßge gives a value of 2.8 h for standard as well as for non-standard conditions, which is in the same order of our TDE deviations caused by variations of the RTML.
The maximum TDE deviation D MAX in CTM and in CMS is observed in the dorsal-ventral axis Y (Fig. 5: SP3, Fig. 6: SP4).
This may be due to the close vicinity of RTML to backbone tissue.A higher value of thermal conductivity of bone tissue of k = 0.75 W/m 2 K compared to the conductivity of other soft tissue types, will result in higher temperature gradient in the line C -nearest backbone parts.Similarly, a higher temperature gradient occurs along the Y: dorsal-ventral SP3-C-SP4 axis in case of CMS due to high conductivity of substrate.Only in CM the axis SP6 -C -SP5 (Z: caudal-cranial) shows the highest D MAX values (see Fig. 7).This may be the case because the CM model lacks the high conductivity substrate the CMS model contains.Other differences between CM and the other models might as well be part of the explanation: The CM shows symmetric anatomical structures with respect to the mid sagittal plane while the CTM does not.Furthermore, the CTM has a much finer discretization than the CM and a varying spatial tissue distribution due to the manual versus the CT-based mesh generation in CM and CTM respectively.Concerning the FE-meshwidth, we found only a minor effect on the TDE difference results as well as [7].Hence, the TDE differences between the CM and the CTM can hardly be explained from differences in the discretization order.
The only difference between CMS and CM is that in CMS, the coarse-meshed human model lies on a discretized structure with wet soil's thermal properties whereas CM consists of a coarse-meshed human model floating CM,CMS compares CM vs. CMS model thus considering the influence of the support structure on the sensitivity of TDE against RTML variation.Like in forensic scenarios heat transfer from the corpse to support with high conductivity and / or heat capacity will increase heat flux through the models support-contact faces.Assuming a support with comparably low thermal conductivity, like a highly insulated mattress, a decrease in heat flux will be the consequence.The influence of different support materials on TDE was investigated in [21,24].Thermal properties of the support were also considered in the Henßge model, where supports with high thermal conductivities are taken into account by increasing the correction factor and vice versa (see e.g.[15]).Our results show that the influence of different supports on TDE sensitivity with respect to RTML variations need to be considered because they may lead to differences in temperature gradients which transform to differences in TDE.
The ambient temperature T A influences the TDE sensitivity against RTML variations in the CMS model in such a way that the TDE deviation increases with decrease in T A .Possibly, this effect is due to the presence of the support structure along the dorsal-ventral axis (Y axis (SP4 -C -SP3)), which increases the rate of heat transfer.The initial temperature of the support structure was always set to the ambient temperature T A .In contrast, looking at the TDE difference D MAX CM,C ™ between D MAX in CM and in CTM, there is nearly no influence of ambient temperature T A .Apparently, the influence of T A on TDE difference is more profound when the support structure is present.
Concerning the influence of RTML variation radius R, the maximum TDE deviation D MAX increases with increase in R in CTM, CMS and CM.This is due to greater Euclidean distance of measurement points from the actual measurement point, which is translated via temperature gradient to greater D MAX .
Figure 10 gives the overview of the maximum TDE deviation by RTML variation evaluated by D MAX for CTM, CMS and CM with respect to ambient temperature T A and measurement variation radius R.
The evaluation of the Qi-located TDE distances D MAX,Qi M between SP R k and C yields the following results: Firstly measurement location errors in the ventral -dorsal axis Y seem to cause the most severe deviations of TDE (see (R1)).
For the coarse model regardless of whether there is a wet soil substrate (M = CMS) or not (M = CM) the D MAX,Qi values essentially decline with rising Q-index I that is for rising cooling times (see (R4)).Since the latter regularity is not repeated by the CT-generated model CTM, we infer this to be an artificial effect caused by the unrealistic high symmetry and regularity of the 'hand-crafted' FE models CMS and CM.Generally the effect of RTML errors on TDE seems to be highest at short times after death (see (R1)).Time evolution of TDE errors due to RTML variations thus correspond to the well known effect, that initial temperature deviations lead to high TDE errors during short times p.m. but are later damped to a low constant basic amount.This is shown in [6] for MHH -TDE but may be generalized.Moreover it is consistent with [7].
There are limitations of the current study.Firstly, we used a simulation-based approach for TDE estimation and RTML variation.A validation of our results on experimental measurements with real corpses was not realized.Measurements using corpses of recently deceased -aside from rare availability of such corpses-are difficult to carry out and especially a reproducible RTML variation is nearly impossible.However, in former studies the CM was calibrated and validated based on some experimental measurements published in [21].Due to the physicsbased approach, both the CM and the CTM should be suitable to carry out parameter variations in terms of RTML variation.Secondly, we only use three different models with normal body mass index.Cooling scenarios with more complex boundary conditions, like adipose corpses, clothing and coverings, etc. were not considered although these factors can influence TDE deviation as CMS depicts.Since varying all those boundary conditions would exceed the scope of one research article, those sensitivity studies are left for future research.

Conclusions
TDE variations caused by RTML deviations may be a source of error in legal medicine TDE since they may reach a magnitude of 2-3 h.A possible consequence for legal medicine applied research could be to investigate methods and devices for reproducible RTML determination.

Fig. 1 a
Fig. 1 a: CT Meshed male corpse model with magnified mesh detail.b: Right half view of CTM model with magnified measurement location showing measurement points C, SPk (k = 1,…,6) in three different octahedra of radii R = 0.5 cm, 1 cm, 2 cm.For the octahedra see Fig. 4.Only the contour lines of the octahedra are shown for clarity

Fig. 2 3
Fig. 2 a: Coarse model lying on a support structure (CMS), b: Coarse model without support structure (CM), c: Detailed view of measurement points C, SP 1 ,….,SP 6 along different radii R = 0.5 cm, 1 cm, 2 cm.For the octahedra see Fig. 4.Only the contour lines of the octahedra are shown for clarity

Fig. 10
Fig. 10 Maximum TDE deviation in D MAX in CTM, CMS and CM model 1 , T 2 ).It is defined in a global version on the whole time interval [t MIN , t MAX ] and in four local ones residing on one of the intersections t Qj ∩ [t MIN , t MAX ] on the time axis each.For all 1 ≤ i ≤ 4 the number of t k lying in t Qj ∩ [t MIN , t MAX ] is denoted by K i.

Table 2
TDE deviation D MAX M in M = CTM at T A = 5° C. Maximum value for fixed R and fixed SP R k: light grey, maximum value for fixed R over all SP R k: fat print, dark grey

Table 3
TDE deviation D MAX M in M = CTM at T A = 15° C. Maximum value for fixed R and fixed SP R k: light grey, maximum value for fixed R over all SP R k: fat print, dark grey

Table 4
TDE deviation D MAX M in M = CTM at T A = 25° C. Maximum value for fixed R and fixed SP R k: light grey, maximum value for fixed R over all SP R k: fat print, dark grey

Table 5
TDE deviation D MAX M in M = CMS at T A = 5° C. Maximum value for fixed R and fixed SP R k: light grey, maximum value for fixed R over all SP R k: fat print, dark grey

Table 6
TDE deviation D MAX M in M = CMS at T A = 15° C. Maximum value for fixed R and fixed SP R k: light grey, maximum value for fixed R over all SP R k: fat print, dark grey

Table 7
TDE deviation D MAX M in M = CMS at T A = 25° C. Maximum value for fixed R and fixed SP R k: light grey, maximum value for fixed R over all SP R k: fat print, dark grey

Table 8
TDE deviation D MAX M in M = CM at T A = 5° C. Maximum value for fixed R and fixed SP R k: light grey, maximum value for fixed R over all SP R k: fat print, dark grey

Table 9
TDE deviation D MAX M in M = CM at T A = 15° C. Maximum value for fixed R and fixed SP R k: light grey, maximum value for fixed R over all SP R k: fat print, dark grey

Table 10
TDE deviation D MAX M in M = CM at T A = 25° C. Maximum value for fixed R and fixed SP R k: light grey, maximum value for fixed R over all SP R k: fat print, dark grey