Double exponential model of corpse cooling under conditions of linearly varying ambient temperature

全文:

ОБОСНОВАНИЕ

Термометрическое определение давности наступления смерти (ДНС) по-прежнему строится с учётом степени охлаждения практически одних лишь только глубоких тканей трупа. Приоритет термометрии ядра тела в диаг-ностике ДНС объясняется более медленным остыванием глубоких тканей, позволяющим увеличить продолжительность посмертного интервала, доступного диагностике, и меньшей подверженностью влиянию различных случайных факторов на процесс охлаждения. Наиболее востребованными диагностическими показателями при этом остаются ректальная и краниоэнцефальная температура, а основой определения ДНС является предложенная в 1962 г. T.K. Marshall и F.E. Hoare феноменологическая модель вида

 $\frac{T-T_a}{T_0-T_a}=\frac{p}{p-k}{e}^{-kt}-\frac{k}{p-k}{e}^{-pt}$, (1)

где Т ― текущая температура ядра тела, °С; Т_а ― температура внешней среды, °С; Т₀ ― начальная температура в момент наступления смерти человека, °С; k ― постоянная охлаждения; p ― постоянная температурного плато; t ― ДНС, ч [1].

Модель (1) охлаждения ядра трупа была построена T.K. Marshall и F.E. Hoare путём введения в закон охлаждения Ньютона-Рихмана [$\frac{dТ}{dt}=-k(T-{Т}_a)$], предназначенный для описания остывания поверхностей физических тел, дополнительной экспоненциальной функции, аппроксимировавшей феномен температурного плато. С учётом этого закон охлаждения ядра трупа при постоянной температуре внешней среды приобрёл вид дифференциального уравнения

$\frac{dТ}{dt}=-k(T-{Т}_a)+С{е}^{-pt}$, (2)

где С ― коэффициент пропорциональности дополнительной экспоненциальной функции, определяемый как

$C=k(T_0-T_a)$. (3)

Решением дифференциального уравнения (2) и явилась феноменологическая модель (1), получившая из-за наличия в своей структуре второй экспоненциальной функции наименование двойной экспоненциальной. Следует отметить, что в строгом смысле двойной экспоненциальной называется функция с наличием второй экспоненты в показателе степени:

$f\left(x\right)=e^{\left(e^x\right)}$.

В этой связи модель охлаждения (1) Marshall-Hoare представляет собой не двойную экспоненциальную, а сумму двух отдельных экспоненциальных функций. Однако термин «двойная экспоненциальная модель» прочно укоренился в судебно-медицинской научной литературе, поэтому при дальнейшем изложении, согласно сложившейся традиции, данным термином будет именоваться сумма двух «простых» экспоненциальных функций.

В 1988 г. немецкий судебный медик C. Henssge предложил собственные методы генерации значений постоянных охлаждения и температурного плато в модели (1) и осуществил масштабную проверку её точности на экспериментальном материале и практических наблюдениях. В результате, модификации C. Henssge позволили определять не только точечные, но и интервальные оценки ДНС с учётом характера одежды на трупе, особенностей его ложа и внешней среды [2, 3].

Наиболее существенный и вследствие этого сильно ограничивающий область применения недостаток двойной экспоненциальной модели (1) характеризуется её корректностью только при постоянной температуре внешней среды. По этой причине закономерными явились многочисленные попытки адаптации модели (1) к условиям изменяющейся внешней температуры. В результате удалось разработать эмпирические модификации двойной экспоненциальной модели, позволившие определять ДНС при однократном дискретном снижении или повышении внешней температуры с последующим сохранением её постоянства [4, 5]. Однако построить феноменологическую модель охлаждения ядра тела для остальных типов возможных изменений внешней температуры, имеющих более актуальное значение, до сих пор так и не удалось.

Вместе с тем основной недостаток двойной экспоненциальной модели (1), связанный с её непригодностью при изменяющейся внешней температуре, может быть преодолен путём обобщения закона (2) на указанные условия охлаждения тела. В целом решать подобные задачи возможно методами прямого или обратного моделирования охлаждения трупа [6]. Однако к уравнению (2) Marshall-Hoare из-за наличия в его структуре дополнительной экспоненциальной функции температурного плато, содержащей в качестве одной из переменных показатель ДНС, применим только метод прямого моделирования. Последний, в свою очередь, также может быть успешно использован только при условии постоянства скорости изменений внешней температуры, т.е. при соответствии изменений температуры внешней среды линейному закону [6]. В связи с тем, что двойная экспоненциальная модель, как и любая её модификация, представляет собой неявно заданную функцию, которую невозможно преобразовать с выведением показателя ДНС в явном виде, актуальным при построении модификаций уравнения (1) является также поиск численного алгоритма их решения с последующей его реализацией в формате компьютерной программы [7].

Цель исследования ― построение математической модели охлаждения ядра трупа на основе феноменологического закона Marshall-Hoare в условиях линейно изменяющейся внешней температуры; поиск численного алгоритма решения построенной модели и последующая его реализация в формате компьютерной программы.

МАТЕРИАЛ И МЕТОДЫ

Методологический дизайн исследования представляет собой математическое моделирование охлаждения ядра трупа в условиях линейно изменяющейся температуры внешней среды, выполненное на базе феноменологичес-кого закона охлаждения Marshall-Hoare и ориентированное на решение задачи определения ДНС по ректальной или краниоэнцефальной температуре.

Математическое моделирование осуществляли прямым аналитическим методом с использованием приложения Wolfram|Alpha. Частные решения построенной математической модели охлаждения ядра трупа получали методом хорд. Вычислительные процедуры, связанные с реализацией итерационного алгоритма метода хорд, производили с использованием приложений Microsoft Excel пакета Office 2016. Визуализацию кривых охлаждения выполняли с помощью приложения Statistica (StatSoft) версии 7.0. Код программы для ЭВМ, реализующей построенную математическую модель охлаждения ядра трупа, а также итерационный алгоритм её решения, составляли на языке программирования C# с использованием приложения Microsoft Visual Studio 2019.

РЕЗУЛЬТАТЫ

Метод построения основанной на законе Marshall-Hoare модели охлаждения ядра тела в условиях линейно изменяющейся внешней температуры связан с заменой в дифференциальном уравнении (2) константы температуры внешней среды на соответствующую линейную функцию и последующим решением полученного дифференциального уравнения.

Пусть труп с начальной температурой в диагностической зоне Т₀ находится в среде с начальной температурой Т_а0, и охлаждение ядра тела протекает согласно уравнению (2) Marshall-Hoare с постоянными коэффициентами k и p. В процессе охлаждения температура внешней среды изменяется линейно по закону

$T_a=T_{a0}+\beta t$, (4)

где β ― скорость изменения внешней температуры, °С/ч.

Необходимо описать динамику температуры ядра тела в указанных условиях.

Для линейно изменяющейся внешней температуры в дифференциальном уравнении (2) следует заменить показатель Т_а, в том числе и в составе константы дополнительной экспоненциальной функции (3), линейным уравнением (4):

$\frac{dТ}{dt}=-k(T-{Т}_{a0}-\beta t)+k\left(T_0-T_{a0}-\beta t\right)e^{-pt}$. (5)

Решением полученного дифференциального уравнения (5) является функция

$\begin{array}{c}T=T_{a0}+\beta \left(t-\frac{1}{k}\right)+\left(T_0-T_{a0}+\frac{\beta }{k}\right)e^{-kt}+\frac{k}{k-p}\times \\ \times \left(\left(T_0-T_{a0}-\frac{\beta \left(kt-pt-1\right)}{k-p}\right)e^{-pt}-\left(T_0-T_{a0}+\frac{\beta }{k-p}\right)e^{-kt}\right).\end{array}$ (6)

Из условия (4) следует, что

$T_{a0}={Т}_{а}-\beta t$. (7)

Подставляя правую часть (7) вместо Т_а0 в формулу (6), окончательно находим феноменологический закон охлаждения ядра тела в условиях линейно изменяющейся температуры внешней среды:

$\begin{array}{c}T\left(t\right)=T_a-\frac{\beta }{k}+\left(T_0-T_a+\beta t+\frac{\beta }{k}\right)e^{-kt}+\frac{k}{k-p}\times \\ \times \left(\left(T_0-T_a+\beta t-\frac{\beta \left(kt-pt-1\right)}{k-p}\right)e^{-pt}-\left(T_0-T_a+\beta t+\frac{\beta }{k-p}\right)e^{-kt}\right).\end{array}$(8)

Классическая модель (1) Marshall-Hoare, найденная при равной нулю скорости изменения внешней температуры, т.е. при её постоянстве, является частным случаем общего уравнения (8). Аналогично закон охлаждения (2) является частным случаем уравнения (5).

Как и в классической модели (1) Marshall-Hoare, неизвестные коэффициенты k и p в уравнении (8) не могут быть однозначно определены путём динамической термометрии трупа. Именно поэтому постоянные охлаждения и температурного плато по-прежнему целесообразно определять в соответствии с рекомендациями С. Henssge [2, 3, 8, 9].

Пример 1

Труп массой 100,0 кг обнаружен в безветренную погоду лежащим на сухой земле в затенённом месте. На трупе два слоя тонкой одежды. Ректальная температура трупа 34,10°С, краниоэнцефальная ― 26,16°С. Температура воздуха на уровне тела на момент его осмотра составила 22,0°С, при этом в предшествующие осмотру 10 ч температура воздуха возрастала по линейному закону со скоростью 0,5°С/ч. Необходимо определить ДНС.

Согласно С. Henssge, корректирующий фактор для трупа массой 70 кг в тонкой двухслойной одежде равен 1,1 [3, 8, 9]. Адаптация корректирующего фактора для ложа трупа не проводится. Адаптацию корректирующего фактора по массе тела находим по формуле $f=\frac{{\left[\frac{-\mathrm{1,2815}}{-\mathrm{3,24596}{e}^{-\mathrm{0,89959}\cdot \mathrm{1,1}}\cdot \left({100}^{-\mathrm{0,625}}-\mathrm{0,028}\right)-\mathrm{0,0354}}\right]}^{\mathrm{1,6}}}{100}=\mathrm{1,06043}$.

Тогда

$k=\mathrm{0,0284}-\mathrm{1,2815}{(\mathrm{1,06043}\cdot 100)}^{-\mathrm{0,625}}=-\mathrm{0,04107}$.

Поскольку модель (8) является обобщением уравнения (1) Marshall-Hoare, то она предполагает положительное значение постоянной охлаждения, поэтому для дальнейших расчетов возьмём последнюю с противоположным знаком.

Так как температура окружающей среды на момент осмотра трупа и в предшествующий ему период была меньше 23,2°С, то константа температурного плато в данном случае равна $p=5k=5\cdot \mathrm{0,04107}=\mathrm{0,205346}$.

Тогда

$\begin{array}{c}0=\left(\mathrm{37,2}-22+\mathrm{0,5}t+\frac{\mathrm{0,5}}{\mathrm{0,041107}}\right)e^{-\mathrm{0,041107}t}+22-\frac{\mathrm{0,5}}{\mathrm{0,041107}}+\frac{\mathrm{0,041107}}{\mathrm{0,041107}-\mathrm{0,205346}}\times \\ \times \left(\begin{array}{l}\left(\mathrm{37,2}-22+\mathrm{0,5}t-\frac{\mathrm{0,5}\left(\mathrm{0,041107}t-\mathrm{0,205346}t-1\right)}{\mathrm{0,041107}-\mathrm{0,205346}}\right)e^{-\mathrm{0,205346}t}-\\ -\left(\mathrm{37,2}-22+\mathrm{0,5}t+\frac{\mathrm{0,041107}}{\mathrm{0,041107}-\mathrm{0,205346}}\right)e^{-\mathrm{0,041107}t}\end{array}\right)-\mathrm{34,1.}\end{array}$ (9)

Решением неявного уравнения (9) является t=9,30 ч. Без учёта изменений внешней температуры оценка ДНС равнялась бы 10,03 ч, а относительная ошибка её определения составила 7,8%.

Пример 2

Найдем ДНС по краниоэнцефальной температуре для данных примера 1. Согласно С. Henssge, константы двойной экспоненциальной модели принимают следующие значения [9]: $k=\mathrm{0,127}$ и $р=\mathrm{1,07}$.

Тогда

$\begin{array}{l}0=22+\left(\mathrm{37,2}-22+\mathrm{0,5}t+\frac{\mathrm{0,5}}{\mathrm{0,127}}\right)e^{-\mathrm{0,127}t}+\frac{\mathrm{0,127}}{\mathrm{0,127}-\mathrm{1,07}}\times \\ \times \left(\begin{array}{l}\left(\mathrm{37,2}-22+\mathrm{0,5}t-\frac{\mathrm{0,5}\left(\mathrm{0,127}t-\mathrm{1,07}t-1\right)}{\mathrm{0,127}-\mathrm{1,07}}\right)e^{-\mathrm{1,07}t}-\\ -\left(\mathrm{37,2}-22+\mathrm{0,5}t+\frac{\mathrm{0,127}}{\mathrm{0,127}-\mathrm{1,07}}\right)e^{-\mathrm{0,127}t}\end{array}\right)-\frac{\mathrm{0,5}}{\mathrm{0,127}}-\mathrm{26,16.}\end{array}$ (10)

Итерационным методом находим, что t=9,30 ч. Без учёта изменений внешней температуры оценка ДНС равнялась бы 11,18 ч, а относительная ошибка её определения составила уже 20,2%.

Геометрическая интерпретация результатов изложенных примеров приведена на рис. 1, на котором температурный тренд внешней среды условно продолжен до 26 ч посмертного периода. Это позволило визуализировать факт соответствия температуры обеих диагностических точек температуре внешней среды после их выравнивания.

Рис. 1. Температурные кривые из примеров 1 и 2.

Примечание. ДНС ― давность наступления смерти.

Fig. 1. Temperature curves from examples 1 and 2.

Note. ДНС ― the postmortem interval, hours.

Уравнение (8) также представляет собой неявно заданную функцию, на основе которой невозможно аналитическое определение ДНС. Однако значение (8) может быть найдено с помощью численных методов. Для этого необходимо преобразовать уравнение (8) путём переноса переменной T в его правую часть и найти численным методом значение ДНС, обращающее полученную функцию вида f (t) = 0 в нуль. Так, подобным образом частные случаи модели (8) в примерах 1 и 2 преобразованы в уравнения (9) и (10).

В качестве численного метода решения двойных экспоненциальных моделей охлаждения трупа был предложен метод касательных Ньютона [7]. Однако при попытках решения методом Ньютона выяснилось, что полученная от (8) неявная функция имеет два корня, один из которых является искомым решением, а второй ― отрицательным числом (рис. 2). При этом метод Ньютона практически при любых начальных значениях t₀ сходится к отрицательному корню уравнения. По этой причине положительный корень уравнений (9), (10) и остальных частных случаев модели (8) было решено находить с помощью иного численного алгоритма ― метода хорд, в рамках которого сущность итерационного процесса определяется формулой $t_{n+1}=t_n-\frac{f\left(t_n\right)}{f\left(t_n\right)-f\left(t_0\right)}\left(t_n-t_0\right)$, $n=\mathrm{0,1,2,}\dots$

где t₀ и t_n ― начальные значения концов отрезка со значениями ДНС, первый из которых с абсциссой t₀ является неподвижным, а второй сходится к положительному корню.

Выбор неподвижного конца хорды определяется формой неявной функции в области положительных значений переменной t. В случаях, когда данная функция на указанном промежутке является выпуклой кверху, в качестве неподвижного следует выбирать правый конец хорды (см. рис. 2). Подобная форма неявных функций имеет место при определении ДНС по ректальной температуре с оценкой констант охлаждения и температурного плато методом Henssge. Если же функция на указанном промежутке является выпуклой книзу, что наблюдается при определении ДНС по краниоэнцефальной температуре, то неподвижным нужно принять левый её конец (рис. 3).

Рис. 2. Геометрия решения уравнения (9) методом хорд для данных примера 1. Показаны первые 3 итерации при t0=10 и t1=0,1 ч.

Fig. 2. The geometry of the solution of equation (9) by the chord method for the data of Example 1. The first 3 iterations are shown at t₀=10 and t₁=0,1 h.

Рис. 3. Геометрия решения уравнения (10) методом хорд для данных примера 2. Показаны первые 5 итераций при t0=0,1 и t1=20 ч.

Fig. 3. The geometry of the solution of equation (10) by the chord method for the data of Example 2. The first 5 iterations are shown at t₀=0,1 and t₁=20 h.

Учитывая, что асимптотическая скорость сходимости метода хорд меньше, чем метода Ньютона, для нахождения корня (9) с точностью до 2 знаков в дробной части результата требуется всего 3, а с точностью до 9 знаков ― 8 итераций.

Метод хорд может быть легко реализован с помощью наиболее распространённого табличного процессора Microsoft Excel. Структура рабочей области в Excel, предназначенной для решения уравнения (8), зависит от метода вычисления констант охлаждения и температурного плато. Общей чертой любой такой рабочей области является наличие в её составе разделов для введения исходных физических величин, вычисления значений констант охлаждения и температурного плато и собственно реализации метода хорд (рис. 4).

Рис. 4. Рабочая область Excel с данными примера 1 для решения функции (9) применительно к ректальной температуре методом хорд. Цветом выделены ячейки, предназначенные для ввода исходных физических величин. Остальные числовые данные генерируются автоматически.

Fig. 4. The Excel workspace with the data of Example 1 for solving the function (9) in relation to rectal temperature by the chord method. The cells intended for the introduction of initial physical quantities are highlighted in color. The remaining numeric data is generated automatically.

Интервальные оценки ДНС в условиях линейно изменяющейся внешней температуры могут быть определены по тем же принципам, что и аналогичные оценки по методу С. Henssge [2, 8, 9]. Применение этого метода оправдано, поскольку все модели С. Henssge являются частными случаями уравнения (8).

Для удобства математическая модель (8) и итерационный алгоритм её решения, а также вычисления интервальных оценок ДНС были реализованы на языке C# в формате программы Warm Bodies MHNH (свидетельство о государственной регистрации № 2021664764).

Значения постоянных охлаждения и температурного плато приложение определяет в соответствии с рекомендациями С. Henssge [8, 9]. Кроме точечных возможно также определение одно- и двусторонних интервальных оценок ДНС при любом уровне доверительной вероятности. Интервальные оценки ДНС приложение определяет путём вычисления остаточного стандартного отклонения на основе регрессионных аппроксимаций [10].

При определении ДНС по ректальной температуре из-за различий значений константы температурного плато и, следовательно, структуры модели (8) в зависимости от внешней температуры в случаях, когда внешняя температура в посмертном периоде сначала была меньше, а затем больше 23,2 °С, или наоборот, приложение Warm Bodies MHNH вычисляет константу температурного плато по формуле p = 7,5k.

Для работы с программой пользователю необходимо выбрать требуемую диагностическую зону трупа, указать её начальную температуру, результаты измерения температуры тела и внешней среды на момент осмотра трупа, величину статистической ошибки и почасовую скорость изменения температуры внешней среды. При выборе ректальной температуры в качестве измеряемого показателя дополнительно требуется указать вес трупа и значение корректирующего фактора для условий охлаждения и ложа трупа массой 70 кг. Адаптация корректирующего фактора по массе трупа осуществляется приложением. При нулевой скорости изменения внешней температуры (β=0) результат работы программы эквивалентен таковому стандартного метода Henssge.

ОБСУЖДЕНИЕ

Основной недостаток феноменологической модели охлаждения Marshall-Hoare связан с её применимостью только в условиях постоянной внешней температуры. Многочисленные исследования, направленные на адаптацию модели (1) к условиям меняющейся температуры окружающей среды, привели лишь к разработке эмпирических модификаций данной модели, предназначенных для определения ДНС при однократном дискретном снижении или повышении внешней температуры с последующим сохранением её постоянства [4, 5]. Подобные изменения условий охлаждения наблюдаются, например, при однократных перемещениях трупа в среду с иной постоянной температурой. Однако каких-либо математических моделей, основанных на феноменологическом законе охлаждения Marshall-Hoare или на предназначенных для глубокой термометрии аппроксимациях закона охлаждения Ньютона-Рихмана и пригодных для использования при монотонных изменениях внешней температуры, предложено не было. Именно поэтому основные усилия исследователей в последние десятилетия сконцентрировались на построении сложных конечно-элементных моделей посмертного теплообмена [11–13]. Предложенные конечно-элементные модели могут учитывать практически все существенные условия охлаждения, в том числе и изменения внешней температуры, однако из-за своей высокой сложности требуют наличия серьёзной физико-математической подготовки и технических навыков, дорогостоящего программного обеспечения и посмертной компьютерной томографии. По этим причинам они не нашли пока широкого применения в экспертной практике.

В рамках настоящего исследования феноменологический закон охлаждения Marshall-Hoare распространён на случай линейных изменений внешней температуры. На основе указанного обобщения построена феноменологическая модель (8) охлаждения ядра трупа в условиях линейно изменяющейся температуры окружающей среды. Поскольку уравнение (8) представляет собой неявно заданную функцию, не допускающую возможность аналитического поиска её корней, предложено находить их численными методами, наиболее подходящим из которых оказался метод хорд. Данный итерационный метод может быть легко воспроизведён с помощью табличного процессора Excel или приложения Warm Bodies MHNH, написанного на языке C#.

Корректность результатов модели (8) сохраняется только при отсутствии таких факторов изменчивости постоянных охлаждения и температурного плато, как инсоляция, естественная или принудительная конвекция, изменения условий охлаждения мёртвого тела. Кроме того, на практике колебания температуры внешней среды вне помещений обычно имеют линейный характер в период времени продолжительностью до 12 ч, реже дольше. Иногда указанные колебания на больших временных промежутках удаётся аппроксимировать линейной регрессией, однако в большинстве случаев модель (8) может быть востребована лишь в первые 12 ч, предшествующие осмотру трупа на месте его обнаружения.

Не исключено, что использование в модели (8) эмпирических коэффициентов, оценённых C. Henssge в условиях постоянной внешней температуры, может повлиять на результаты определения ДНС при значительных изменениях температуры окружающей среды. Обнаружение подобного влияния потребует адаптации процедуры задания коэффициентов модели (8) для той скорости изменений внешней температуры и их перепадов, при которых оно будет выявлено.

ЗАКЛЮЧЕНИЕ

Разработана математическая модель охлаждения ядра трупа в условиях линейно изменяющейся температуры внешней среды, основанная на феноменологическом законе охлаждения Marshall-Hoare. Классическая двойная экспоненциальная модель Marshall-Hoare является частным случаем новой модели, реализуемым при нулевой скорости изменения внешней температуры. В качестве численного алгоритма решения разработанной математической модели предложен метод хорд.

Разработанная математическая модель и итерационный алгоритм её решения, а также процедуры вычисления интервальных оценок ДНС реализованы на языке C# в формате компьютерной программы Warm Bodies MHNH.

Предложенный метод и реализующую его компьютерную программу целесообразно использовать в судебно-медицинской экспертной практике при определении ДНС по ректальной или краниоэнцефальной температуре трупа в условиях линейно изменяющейся температуры внешней среды.

ДОПОЛНИТЕЛЬНО

Источник финансирования. Автор заявляет об отсутствии внешнего финансирования при проведении исследования.

Конфликт интересов. Автор декларирует отсутствие явных и потенциальных конфликтов интересов, связанных с публикацией настоящей статьи.

Вклад авторов. Автор подтверждает соответствие своего авторства международным критериям ICMJE (разработка концепции, проведение исследования и подготовка статьи, одобрение финальной версии перед публикацией).

ADDITIONAL INFORMATION

Funding source. This study was not supported by any external sources of funding.

Competing interests. The author declare that they have no competing interests.

Authors’ contribution. The author made a substantial contribution to the conception of the work, acquisition, analysis, interpretation of data for the work, drafting and revising the work, final approval of the version to be published and agree to be accountable for all aspects of the work.

作者简介

German Nedugov

编辑信件的主要联系方式.
Email: nedugovh@mail.ru
ORCID iD: 0000-0002-7380-3766
SPIN 代码: 3828-8091
Scopus 作者 ID: 25947646500
Researcher ID: ABH-5590-2020

MD, Dr. Sci. (Med.), Associate Professor

参考

补充文件