Compensation for Aberrations When Focusing Ultrasound Through the Skull Based on CT and MRI Data

Толық мәтін

ВВЕДЕНИЕ

В современной медицине наиболее предпочтительными становятся малоинвазивные и неинвазивные методы хирургии. Одним из таких методов является ультразвуковая нейрохирургия. Метод заключается в фокусировке мощного ультразвука в заданные участки мозга через невскрытый череп, что приводит к локальному нагреву и разрушению тканей [1–3]. Осуществление фокусировки ультразвука в мозг через кости черепа стало возможным благодаря развитию методов визуализации, таких как компьютерная и магнитно-резонансная томография (КТ и МРТ), совершенствованию вычислительных методов и появлению многоэлементных фазированных ультразвуковых решеток [4, 5].

Из-за неоднородности внутренней структуры, толщины и геометрии поверхностей черепа, а также сильного отличия его акустических свойств от прилегающих сред (кожи и мозга) происходит искажение фокусировки ультразвукового пучка, поэтому на практике для проведения операций необходимо компенсировать возникающие аберрации путем подбора соответствующих фаз на элементах терапевтической решетки [5]. В настоящее время операции с использованием мощного фокусированного ультразвука проводятся с применением клинических систем ExAblate (InSightec Ltd., Tirat Carmel, Israel) в виде полусферической фазированной решетки в два этапа. На первом этапе делается рентгеновская компьютерная томография головы пациента. Контраст изображения КТ характеризуется значениями по шкале Хаунсфилда, из которых могут быть рассчитаны плотность и скорость звука в каждой точке изображения [6, 7]. Данные КТ далее используются для первичного отбора пациентов и для расчета компенсирующих фаз на элементах решетки. Отбор пациентов происходит в основном исходя из параметра SDR (Skull Density Ratio), который характеризует степень однородности внутренней структуры костей черепа и определяется как отношение средней плотности губчатого вещества к средней плотности компактного вещества в единицах Хаунсфилда [8]. Для пациентов с высоким показателем SDR (SDR ≥ 0.45) облучение наиболее эффективно и безопасно. На втором этапе для проведения операции пациента помещают в аппарат МРТ с расположенным в нем ультразвуковым излучателем. Перед началом ультразвукового облучения предварительно полученное трехмерное КТ-изображение головы пациента совмещается в режиме реального времени с МРТ-изображением, и дальнейшая операция проводится с использованием МРТ-визуализации для контроля температуры в области воздействия [5].

Необходимые процедуры предварительного получения данных КТ, а также совмещения КТ- и МРТ-изображений перед ультразвуковым облучением увеличивают время подготовки к операции и усложняют ее проведение. Возможность первичного отбора пациентов на основе данных МРТ может позволить избежать воздействия ионизирующего излучения на пациентов, а проведение операции под МРТ-контролем с компенсацией аберраций в режиме реального времени – избежать трудностей с совмещением данных КТ и МРТ. Поэтому одной из практически важных задач современной ультразвуковой нейрохирургии является исследование возможности проведения компенсации аберраций по данным МРТ. Основная сложность заключается в том, что в отличие от КТ, МРТ-изображения не позволяют напрямую получить количественные данные об акустических параметрах тканей головы, а именно плотности тканей и скорости звука в них.

Можно выделить два основных подхода, которые предпринимались для решения данной задачи. Первый подход использует возможности МРТ, основанные на регистрации ультракороткого или нулевого времени эха T2 компонент [9–11]. Развитие данных методов регистрации позволяет визуализировать внутреннюю структуру кости черепа и найти корреляцию между интенсивностью МРТ-изображения и шкалой Хаунсфилда [10, 12]. Модели, построенные на основе таких данных МРТ, учитывают неоднородную структуру черепа, при этом мягким тканям (кожа и мозг) присваивается постоянное одинаковое значение по шкале Хаунсфилда [9, 10] либо рассматривается модель черепа в воде [11]. Другим подходом является использование машинного обучения и нейросетей для построения КТ-изображения в единицах Хаунсфилда по интенсивности МРТ-изображения [13–16]. В данном подходе для обучения нейросетей используются T1-взвешенные МРТ-изображения и МРТ-изображения с ультракоротким и нулевым временем эха в комплекте с КТ-изображением того же пациента [16].

На сегодняшний день исследование возможности компенсации аберраций с использованием акустической модели, построенной на основе данных МРТ, проводилось только с использованием лучевого метода для полусферической решетки системы ExAblate и моделей черепа в воде. Наблюдались отличия в результате проведения компенсации аберраций по моделям МРТ и КТ, однако какие именно параметры моделей и соответствующие волновые эффекты определяли эти отличия, не анализировалось [9]. Таким образом, отдельным важным вопросом является исследование влияния различных характеристик моделей, таких как геометрические отличия восстановленных поверхностей черепа и учет внутренних неоднородностей различных костей черепа, на искажение фокусированного пучка и результат коррекции.

В данной работе рассмотрены данные обследования реального пациента, которому были проведены КТ- и МРТ-сканирования головного мозга. В качестве излучателя использовалась модель 256-элементной решетки с рабочей частотой 1 МГц и абсолютно плотным мозаичным заполнением поверхности элементами [17]. Коррекция аберраций проводилась путем численного моделирования распространения волны от виртуального точечного источника, расположенного в целевой точке, по направлению к поверхности решетки, – т. н. “дифракционный” метод коррекции [18]. Основной целью работы являлось сравнительное исследование возможностей проведения компенсации аберраций по моделям, построенным на основе данных КТ и МРТ. Кроме того, проанализировано влияние различий во внутренней структуре и геометрии различных тканей в построенных моделях на результат компенсации аберраций. Наконец, практической целью работы являлось составление рекомендаций по использованию моделей на основе КТ и МРТ для планирования ультразвукового облучения и проведения коррекции во время операции.

ТЕОРЕТИЧЕСКАЯ МОДЕЛЬ

Построение трехмерных акустических моделей головы человека по данным компьютерной и магнитно-резонансной томографии

Для построения трехмерных акустических моделей головы пациента использовались анонимизированные данные КТ и МРТ одного пациента, предоставленные Медицинским научно-образовательным центром МГУ имени М.В. Ломоносова. Трехмерное КТ-изображение включало в себя набор из 586 аксиальных срезов с разрешением 0.4 × 0.4 × 0.3 мм. Интенсивность КТ-изображений характеризуется яркостью пикселей в единицах Хаунсфилда (HU), которые показывают, насколько рентгеновское излучение ослабляется при прохождении через ткань. В частности, значение HU (Hounsfield Units) для воды равно нулю. Было показано, что единицы HU определяются плотностью ткани, а плотность коррелирует со скоростью звука в ней [6, 7]. Таким образом, из данных КТ можно получить информацию об акустических свойствах ткани (плотности и скорости продольных акустических волн) в каждом пикселе изображения.

Чтобы задать коэффициент поглощения, необходимо сегментировать КТ-изображения и выделить основные сегменты – кожу, череп и мозг, коэффициент поглощения для которых известен из справочной литературы [19]. Сегментация КТ-изображений в данной работе была проведена с помощью открытого программного обеспечения Slicer (slicer.org) [20]. Благодаря сильно различающимся значениям пикселей по шкале Хаунсфилда в кости и мягких тканях, к КТ-изображению возможно было применить метод пороговой сегментации. С помощью этого метода и дополнительных функций программы Slicer (flood filling, islands, logic operations) в отдельные сегменты были выделены мозг, череп, кожа и согласующая среда (вода) (рис. 1).

 

Рис. 1. Аксиальный срез компьютерной томографии головы человека (а) – до проведения сегментации, (б) – до проведения сегментации с другим контрастом, (в) – после сегментации в программе Slicer. Сегменты отмечены цифрами: 1 – внешняя среда (вода), 2 – кожа, 3 – череп, 4 – мозг.

 

Данные МРТ представляли собой T1-взвешенное трехмерное изображение с набором из 176 аксиальных срезов с разрешением 1 × 1 × 1 мм. В отличие от КТ-изображения, МРТ предоставляет только качественную визуальную информацию о тканях и органах. Сегментация МРТ-изображения была проведена в два этапа. На первом этапе изображение сегментировалось в программе открытого доступа BrainSuit21a (brainsuite.org), однако из-за ограниченного разрешения МРТ-изображения и анатомических особенностей строения головы человека (наличия зрительных нервов и т.д.) выделение границ черепа с использованием детектора Марра–Хилдретта, используемого в программе BrainSuit21a, дает неточности [21, 22]. Поэтому, для более точного определения границ тканей, сегментация была доработана вручную в программе Slicer, ориентируясь на визуальные границы анатомических структур (рис. 2) [20].

 

Рис. 2. Аксиальный срез магнитно-резонансной томографии головы человека (а) – до сегментации; (б) – после сегментации в программе BrainSuit; (в) – после доработки в программе Slicer. Цвета сегментов соответствуют таковым на сегментации КТ.

 

Полученные в формате nrrd сегментации далее загружались в программу MATLAB (mathworks.com) и интерполировались на расчетные сетки, одинаковые для КТ- и МРТ-изображений. Из-за различного положения пациента при проведении магнитно-резонансной и компьютерной томографий перед интерполяцией МРТ-изображение было дополнительно повернуто на 13.5° относительно оси Ox и на 1.25° относительно оси Oy, чтобы обеспечить наилучшее совмещение поверхностей черепов в обеих моделях (рис. 3). Ноль расчетной сетки по оси z располагался на макушке (наивысшей точке черепа), нули по осям x и y были выбраны в середине осей, симметрично относительно макушки. Шаг расчетной сетки составлял Δx = Δy = Δz = 0.5 мм. Для данных КТ интерполяция на расчетную сетку проводилась с помощью функции griddedInterpolant программы MATLAB. Для данных МРТ использовалась функция scatteredInterpolant и, поскольку размер вокселя превышал в два раза шаг сетки, каждый воксель делился пополам в каждом из трех направлений, при этом тип среды каждого вокселя оставался прежним. В результате были получены две трехмерные матрицы, размером 236 × 236 × 180 точек, представляющие сегментации КТ- и МРТ-изображений, в каждой ячейке которых хранился номер соответствующего сегмента (1 – “вода”, 2 – “кожа”, 3 – “череп”, 4 – “мозг”). На основе КТ-изображения была также создана матрица 236 × 236 × 180 точек, содержащая в каждой ячейке значение по шкале Хаунсфилда. Таким образом, были получены три акустические модели головы человека.

 

Рис. 3. Результат совмещения трехмерных моделей головы человека (слой черепа), построенных на основе данных КТ (синий) и данных МРТ (красный); (а) – внешние поверхности, (б) – внутренние поверхности.

 

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

Для оценки степени однородности черепа были посчитаны два параметра: SS (Skull Score) и SDR (Skull Density Ratio) [23]. Для расчета параметра SS были проведены лучи из геометрического центра каждого элемента к точке фокуса. На каждом луче в области, соответствующей сегменту черепа, выбирались минимальное и максимальное значения по шкале Хаунсфилда и вычислялось их отношение. Параметр SS вычислялся как среднее полученных отношений по всем лучам:

$SS=\frac{1}{N}\sum _{i=1}^N\frac{H{U}_{\mathrm{мин}i}}{H{U}_{\text{макс}i}},$ (1)

где N – количество элементов излучателя, i – номер элемента, HU_минi – минимальное значение по шкале Хаунсфилда на i-ом луче, HU_максi – максимальное значение по шкале Хаунсфилда на i-ом луче.

Кроме этого, был также вычислен параметр SDR, как отношение среднего значения единиц Хаунсфилда вокселей, соответствующих трабекулярной (губчатой) кости, к среднему значению единиц Хаунсфилда вокселей, соответствующих кортикальной кости, в конусе облучения, образованном лучами, исходящими от краев излучателя в точку фокуса:

$SDR=\frac{H{U}_{\text{т}}}{H{U}_{\text{к}}},$ (2)

где HU_т – среднее значение по шкале Хаунсфилда для губчатой кости, HU_к – среднее значение по шкале Хаунсфилда для кортикальной кости. Границы между кортикальной и губчатой костью в единицах Хаунсфилда определялись с помощью алгоритма Оцу [23, 24].

Методы расчета поля и компенсации аберраций

В работе рассматривалась компьютерная модель 256-элементной решетки, предложенная в Лаборатории медицинского и промышленного ультразвука МГУ [17]. Решетка имеет форму сферического сегмента с фокусным расстоянием F = 200 мм, апертурой D = 200 мм и рабочей частотой 1 МГц. Элементы решетки представляют собой многоугольники равной площади, расположенные на поверхности излучателя случайным образом, при этом обеспечивается абсолютно плотное заполнение поверхности элементами (рис. 4а) [17].

 

Рис. 4. (а) – Схема рассматриваемой в работе 256-элементной решетки; диаметр D = 200 см, радиус кривизны F = 200 см, рабочая частота f = 1 МГц, зазор между элементами 0.5 мм. Схемы (б) – расчета поля при облучении мозга и (в) – компенсации аберраций. Стрелками показан порядок проведения расчетов: с использованием интеграла Рэлея и псевдоспектрального метода моделирования распространения акустических волн в неоднородных средах. Фокусировка проводилась в центр мозга.

 

Численный алгоритм моделирования фокусировки ультразвука в мозг состоял из комбинации двух методов. На первом этапе акустическое поле излучателя рассчитывалось с помощью интеграла Рэлея в однородной среде (воде) от поверхности решетки до горизонтальной плоскости (xy), расположенной в воде вблизи поверхности черепа (рис. 4б) [25]. Для ускорения расчетов в условиях данной задачи для расчета интеграла Рэлея использовалась суперпозиция аналитических решений для каждого из элементов решетки следующим образом [18, 25, 26]. Многоугольные элементы излучателя предварительно разбивались на прямоугольные треугольники, для каждого из которых можно выписать аналитическое решение в дальнем поле [25]. Таким образом, поле акустического давления на плоскости xy, перпендикулярной оси решетки и расположенной вблизи поверхности черепа, рассчитывалось как сумма полей подэлементов в форме прямоугольных треугольников, на которые были разбиты элементы излучателя.

Далее для расчета поля внутри головы человека использовался псевдоспектральный метод моделирования распространения акустических волн в неоднородной среде, реализованный с помощью программного пакета k-Wave (k-wave.org) [18, 27–29]. Как указано выше, для моделирования использовалась расчетная сетка с пространственными шагами Dx = Dy = Dz = 0.5 мм. Шаг по времени расчетной сетки выбирался в соответствии с критерием Куранта–Фридрихса–Леви и удовлетворял условию $CFL\equiv c_{\max }\Delta t/\Delta x$, где CFL = 0.1, c_max – максимальная скорость звука в рассматриваемой модели. Полученное на первом этапе распределение акустического давления на плоскости вблизи головы использовалось в качестве граничного условия. Расчетная модель k-Wave учитывает такие эффекты, как наличие неоднородностей акустических параметров, рефракция, отражение и поглощение акустических волн. Эффект генерации сдвиговых волн в черепе не учитывался, так как ранее было показано, что при падении акустического пучка на поверхность черепа, близком к нормальному, генерация сдвиговых волн слабо влияет на искажение ультразвукового пучка [26]. В то же время пренебрежение данным эффектом существенно облегчает и ускоряет расчет в неоднородной области головы пациента.

Расчет поля без компенсации аберраций проводился путем облучения головы пациента с одинаковыми фазами и амплитудами на элементах излучателя. Алгоритм компенсации аберраций также состоял из комбинации двух этапов, в которых использовались различные численные методы (рис. 4в). На первом этапе моделировалось распространение сферической волны из точки фокуса до плоскости xy. Моделирование проводилось с использованием программного комплекса k-Wave [18, 26–29]. На втором этапе решалась обратная задача выбора набора фаз на элементах излучателя, который наилучшим образом воспроизводил бы поле, распространенное обратно на плоскость. Результат расчета поля на плоскости xy был представлен в виде дискретной матрицы с известными комплексными амплитудами давления в узлах сетки. Далее матрица была преобразована в вектор комплексных амплитуд давления на плоскости xy. Каждый элемент данного вектора представлялся как комбинация аналитических решений для дальнего поля элементов излучателя с варьируемой фазой. В дальнейшем полученная система решалась методом наименьших квадратов [18]. В результате были получены корректирующие аберрации фазы для каждого элемента 256-элементного излучателя. Полученные фазы инвертировались, амплитуда устанавливалась одинаковой и равной единице, и далее расчет поля проводился уже с компенсацией аберраций.

Для оценки степени искажения пучка и качества фокусировки были использованы следующие параметры: максимум амплитуды давления, нормированной на амплитуду давления на элементе излучателя p_A/p₀, его положение относительно центра излучателя и эффективность фокусировки [11, 16, 30–33]. Эффективность оценивалась как отношение интегралов интенсивности полученного распределения к распределению поля в воде, взятых в фокальном объеме, определяемом по уровню –8 дБ для случая фокусировки в воде [34]:

$\epsilon =\frac{{\displaystyle \underset{W}{\int }I\left(\mathbf{r}\right)dV}}{{\displaystyle \underset{W}{\int }{I}_{\text{вода}}\left(\mathbf{r}\right)dV}}.$ (3)

Здесь I(r) – интенсивность ультразвука, W – фокальный объем, определяемый по уровню –8 дБ для случая фокусировки в воде.

РЕЗУЛЬТАТЫ

Акустические модели головы человека, построенные по данным компьютерной и магнитно-резонансной томографии

Далее приводятся результаты, полученные для трех основных акустических моделей, построенных по данным КТ и МРТ. Кроме того, для исследования влияния неоднородностей тканей на возможности компенсации аберраций при фокусировке ультразвукового пучка в мозг дополнительно рассматривались несколько модификаций акустической модели, построенной на основе данных КТ.

В первую очередь была проведена оценка точности совмещения изображений КТ и МРТ. Для этого сравнивались толщины и разность толщин черепов в моделях КТ и МРТ, рассчитанные по лучам. Средняя толщина черепа в модели КТ составила 7.33 ± 1.03 мм, в модели МРТ – 7.28 ± 1.34 мм. Отличие толщин, рассчитанных по пути лучей от фокуса к центрам элементов решетки в сегменте череп, в моделях КТ и МРТ лежало в интервале от 0.14 до 1.76 мм. Были найдены также параметры SS и SDR; для данного черепа и данной геометрии облучения оба параметра были одинаковы SS = SDR = 0.51. На практике, черепа со значениями SDR (или SS) выше 0.45 считаются пригодными для проведения транскраниального облучения ультразвуком. Для таких черепов внутренняя структура кости черепа является достаточно однородной и не вносит существенного затухания и сильных искажений в распространение ультразвукового пучка [8].

Для разделения эффектов аберраций, связанных с геометрическими и внутренними неоднородностями черепа, были построены следующие три основные модели головы пациента. Первая модель (неоднородная модель КТ, нКТ) учитывала неоднородности внутренней структуры черепа (рис. 5а). Из единиц Хаунсфилда были рассчитаны плотность и скорость звука для каждого вокселя изображения, а коэффициенты поглощения задавались постоянными для каждой ткани на основе сегментации и представлены в табл. 1 [6, 7, 19]. Вторая модель (однородная модель КТ, оКТ) также была построена на основе данных КТ, но имела постоянные значения акустических параметров для каждого типа ткани (рис. 5б). Данные значения рассчитаны как средние внутри конуса облучения. В третьей модели (однородная модель МРТ) (рис. 5в) каждому сегменту были присвоены постоянные значения, рассчитанные как средние внутри конуса облучения из неоднородной модели КТ. Таким образом, модели оКТ и МРТ отличались лишь геометрическими параметрами черепа. Значения акустических параметров для однородных моделей КТ и МРТ представлены в табл. 2. Во всех моделях внешней контактной средой между головой и излучателем считалась вода. Плотность и скорость звука в воде выбирались равными ρ₀ = 1000 кг/м³ и c₀ = 1500 м/с. Коэффициент поглощения для воды полагался равным нулю. Заметим, что вторая модель оКТ по сути является переходной между нКТ и МРТ: она не учитывает внутренние неоднородности сегментов подобно МРТ, однако геометрически не отличается от нКТ. Такое построение модели позволяет отделить эффект внутренних неоднородностей тканей и исследовать его отдельно.

 

Рис. 5. (а) – Распределение плотности в кортикальной плоскости zy в неоднородной модели КТ. Распределения сегментов на примере кортикального среза zy (б) – в однородной модели КТ и (в) – в модели МРТ, цифрами обозначены индексы каждого сегмента.

 

Таблица 1. Диапазоны значений плотности ρ₀ и скорости продольных волн c_l, рассчитанные из единиц Хаунсфилда КТ-изображения, и значения коэффициента поглощения α на частоте 1 МГц, взятые из справочной литературы, для каждого типа тканей [6, 7, 19]

	ρ0, кг/м3	cl, м/с	α, дБ/см
Кожа	430–1240	1090–1780	1.84
Череп	1060–2500	1575–3100	8.83
Мозг	490–1300	940–1800	0.21

 

Таблица 2. Значения плотности ρ₀, скорости звука продольных волн c_l и коэффициента поглощения α на частоте 1 МГц в однородных моделях КТ и МРТ [19]

	ρ0, кг/м3	cl, м/с	α, дБ/см
Кожа	929	1431	1.84
Череп	1732	2331	8.83
Мозг	1058	1575	0.21

 

Для анализа влияния неоднородностей тканей дополнительно были построены еще три модификации неоднородной модели КТ. В каждой из этих моделей одна из сред (кожа, череп, мозг) считалась однородной, и ей присваивались постоянные значения плотности и скорости звука из табл. 2. Остальные два сегмента сохраняли свою неоднородную структуру, и их акустические параметры рассчитывались из единиц Хаунсфилда.

Искажение акустического поля и компенсация аберраций в рамках одной модели

Для исследования искажения ультразвукового пучка при фокусировке в мозг были проведены расчеты поля при фокусировке в центр мозга на глубину 60 мм, считая от внутренней поверхности черепа, для каждой из трех основных моделей: однородной модели МРТ, однородной модели КТ (оКТ) и неоднородной модели КТ (нКТ). Были проанализированы коэффициент усиления по давлению в максимуме поля относительно характерного давления на поверхности излучателя, p_max/p₀, его положение относительно центра кривизны (геометрического фокуса) излучателя и эффективность фокусировки в объеме распределения.

На рис. 6 представлены распределения нормированной амплитуды давления p_A/p₀ в кортикальной плоскости zy для трех акустических моделей: неоднородной модели КТ (нКТ), однородной модели КТ (оКТ) и модели МРТ. Коэффициент усиления в максимуме поля в модели нКТ составил p_A/p₀ = 29.06, эффективность ε = 8.9%. Смещение фокуса в данной модели было наименьшим и составляет –0.5 мм вдоль оси Ox и –1.5 мм вдоль оси Oz, вдоль оси Oy смещение отсутствует. Для модели оКТ соответствующие значения равны p_A/p₀ = 27.15 и ε = 5.7%, смещение в поперечном направлении совпадает со случаем модели нКТ и составляет –0.5 мм вдоль оси Ox и отсутствует вдоль оси Oy. Вдоль оси Oz смещение больше, чем в модели нКТ, и составляет –2.5 мм. Для модели МРТ: p_A/p₀ = 30.85 и ε = 5.1%. В модели МРТ увеличивается смещение вдоль поперечных осей и составляет –1.5 мм вдоль оси Ox и –0.5 мм вдоль оси Oy, вдоль оси Oz смещение составляет –2.0 мм. Форма пучка после прохождения через череп также различна для разных моделей. Однако при этом амплитудные значения максимума поля отличаются не более чем на 10%. Различие эффективностей обусловлено в первую очередь различным смещением максимума поля относительно геометрического центра кривизны излучателя: эффективность уменьшается при увеличении смещения фокальной области.

 

Рис. 6. Распределения нормированной амплитуды давления p_A/p₀ в плоскости zy в случае облучения центра мозга без компенсации аберраций для (а) – неоднородной модели КТ; (б) – однородной модели КТ; (в) – однородной модели МРТ. Обозначение модели, для которой проводился расчет поля: нКТ – неоднородная модель КТ, оКТ – однородная модель КТ, МРТ – модель МРТ.

 

Компенсация аберраций далее рассчитывалась по той же модели, которая в дальнейшем использовалась при облучении. На рис. 7 представлены распределения нормированной амплитуды давления p_A/p₀ в плоскости zy для всех трех моделей. Надписи вида “нКТ/нКТ” в верхнем правом углу рисунков соответствуют моделям, для которых были проведены коррекция аберраций и расчет поля. По сравнению со случаем облучения без компенсации аберраций (рис. 6) видно, что во всех случаях коррекция позволяет добиться точной фокусировки, избежать образования побочных максимумов и “размытия” фокуса.

 

Рис. 7. Распределения нормированной амплитуды давления p_A/p₀ в плоскости zy в случае облучения центра мозга после компенсации аберраций для (а) – неоднородной модели КТ; (б) – однородной модели КТ; (в) – однородной модели МРТ. Обозначения моделей, используемых для компенсации аберраций и расчета поля: нКТ – неоднородная модель КТ, оКТ – однородная модель КТ, МРТ – модель МРТ.

 

Как следует из результатов, приведенных на рис. 7, после компенсации аберраций во всех моделях КТ положение максимума амплитуды поля p_max/p₀ совпадает с геометрическим центром кривизны излучателя. Для модели МРТ отсутствует сдвиг фокуса относительно поперечных осей, а сдвиг вдоль оси излучателя z составляет –0.5 мм, что порядка размера одного вокселя. Компенсация аберраций позволяет значительно повысить амплитуду давления в фокусе и эффективность, а также получить узкую фокальную перетяжку для всех моделей, шириной не более 2 мм по уровню –3 дБ от максимума, как и при фокусировке в воде. После коррекций аберраций коэффициент усиления в фокусе и эффективность составляют p_F/p₀ = 40.33, ε = 17.5% для модели нКТ, p_F/p₀ = 41.81, ε = 17.0% для модели оКТ и p_F/p₀ = 43.43, ε = 18.8% для модели МРТ. Различие значения давления в фокусе и эффективности для данных моделей также не превосходит 10%. Таким образом, модель МРТ может быть использована для прогнозирования амплитуды и структуры акустического поля в фокальной области пучка при планировании облучения.

Компенсация аберраций по акустическим моделям с однородными тканями

Результаты, приведенные в этом разделе, соответствуют имитации проведения реальной операции, когда расчет поля происходит с использованием наиболее точной модели нКТ, однако коррекция при этом осуществляется с использованием упрощенных моделей оКТ и МРТ. Были проведены расчеты компенсации аберраций с помощью расчета фаз на элементах излучателя по моделям оКТ и МРТ, в дальнейшем фазы инвертировались и при облучении использовалась модель нКТ, наиболее близкая к реальному случаю. Модели оКТ и нКТ отличаются внутренней структурой тканей: модель оКТ не учитывает внутренние неоднородности тканей. Расчет компенсации по модели оКТ и дальнейшее облучение с фазами по модели нКТ позволяет определить влияние внутренних неоднородностей на возможности компенсации аберраций. Модель МРТ, кроме однородной внутренней структуры тканей, имеет также геометрические отличия от модели нКТ. Использование фаз, рассчитанных по модели МРТ, для компенсации аберраций позволяет также оценить влияние погрешности в определении границ тканей на компенсацию аберраций.

Результаты данных расчетов представлены на рис. 8б и сравниваются со случаем проведения компенсации аберраций по неоднородной модели КТ (рис. 8а). Видно, что компенсация аберраций, проведенная по однородной модели КТ, позволяет достичь максимума давления p_A/p₀ = 34.79, что составляет 87% от давления в фокусе при компенсации по неоднородной модели КТ. Эффективность при этом составляет ε = 12.6%, что ниже, чем при компенсации по неоднородной модели КТ (ε = 17.5%). Кроме того, появляется смещение максимума поля относительно геометрического центра кривизны излучателя на 2.0 мм вдоль оси Oz (рис. 8б). В случае компенсации аберрации по модели МРТ, максимум давления равен p_A/p₀ = 26.56 и составляет 65% от случая компенсации по неоднородной модели КТ. Эффективность при компенсации по модели МРТ снижается до 5.3% и возникает смещение вдоль оси Ox, смещение вдоль других осей отсутствует (рис. 8в). Таким образом, для получения эффективной и точной фокусировки необходимо учитывать как внутреннюю структуру тканей, так и наиболее точно восстанавливать их геометрию. При этом неоднородная модель КТ и однородная модель КТ имеют одинаковую геометрию тканей и различаются только их внутренней структурой. Расчет компенсации аберраций по модели, не учитывающей неоднородности тканей, приводит к смещению максимума поля и снижению амплитуды давления в максимуме и эффективности. Модель МРТ и однородная модель КТ имеют одинаковые акустические параметры тканей, но из-за ошибок в сегментации МРТ-изображения отличаются геометрией тканей. Ошибка в учете геометрии также приводит к смещению фокуса и к еще большему снижению в нем амплитуды давления и эффективности.

 

Рис. 8. Распределения нормированной амплитуды давления p_A/p₀ в плоскости zy в случае облучения центра мозга после компенсации аберраций с фазами рассчитанными по (а) – неоднородной модели КТ; (б) – однородной модели КТ; (в) – однородной модели МРТ. Обозначения моделей, используемых для компенсации аберраций и расчета поля: нКТ – неоднородная модель КТ, оКТ – однородная модель КТ, МРТ – модель МРТ.

 

Наконец, было проведено исследование влияния неоднородностей различных тканей на возможности компенсации аберраций. Для этого коррекция аберраций проводилась по моделям, в каждой из которых одна из тканей считалась однородной, а остальные неоднородными. Всего было рассмотрено три модели: с однородной кожей и неоднородными черепом и мозгом, с однородным черепом и неоднородными кожей и мозгом, с однородным мозгом и неоднородными кожей и черепом. Наибольшее смещение максимума поля относительно геометрического центра кривизны излучателя характерно для модели с однородным черепом и составляет 2 мм вдоль оси излучателя Oz. Для модели с однородным черепом характерна также наименьшая амплитуда давления в фокусе p_A/p₀ = 37.11 и наименьшая эффективность ε = 13.6%. Таким образом, проведение компенсации аберраций по модели, не учитывающей неоднородную структуру ткани, приводит к смещению максимума поля относительно целевой точки и уменьшению максимума давления и эффективности. При этом наибольшее влияние вносят неоднородности черепа.

ЗАКЛЮЧЕНИЕ

В данной работе, на примере построенных акустических моделей головы конкретного пациента, было проведено сравнение возможностей использования данных МРТ вместо данных КТ для планирования транскраниального облучения глубоких структур головного мозга с использованием мощного фокусированного ультразвука. Было показано, что при планировании облучения различия между тремя различными акустическими моделями (нКТ, оКТ, МРТ) не являлись существенными. Так, при моделировании облучения без компенсации аберраций различие основных оцениваемых параметров акустического поля, таких как амплитуда в максимуме давления и эффективность, не превышало 10%. Однако наблюдались различие в форме искаженной фокальной перетяжки и небольшое, около 1 мм, смещение точки максимума для разных моделей.

Алгоритм компенсации аберраций, реализованный в рамках одной модели (расчет корректирующих фаз и последующий расчет поля с использованием одной и той же модели головы), позволяет добиться точной фокусировки ультразвука, узкой фокальной перетяжки и повышения давления в фокусе на 40–54% во всех моделях. Полученные для разных моделей результаты также отличаются незначительно. Амплитудные значения давления в фокусе после коррекции и эффективность отличаются не более, чем на 10% при использовании следующих комбинаций “коррекция/последующий расчет” нКТ/нКТ, оКТ/оКТ, МРТ/МРТ. Таким образом, любая из трех моделей может быть использована отдельно для прогнозирования результата фокусировки ультразвука в мозг, предсказания расположения максимума поля и амплитуды давления в фокусе. Отметим, что данный результат применим к черепам с высоким значением SDR (SDR ≥ 0.45), которые обладают плотной структурой губчатой кости черепа. При этом в рассмотренной задаче используется излучатель, конус облучения которого попадает на область макушки и не затрагивает наиболее неоднородные участки костей черепа. Таким образом, данная геометрия излучателя позволяет сократить влияние неоднородностей внутренней структуры черепа на искажение и ослабление ультразвукового пучка.

При моделировании применения различных моделей головы человека для коррекции аберраций в реалистичных условиях, акустические модели нКТ, оКТ и МРТ показали значительное различие в качестве коррекции. Для демонстрации были проведены расчеты компенсации аберраций по однородным моделям КТ и МРТ, после чего облучение проводилось по неоднородной модели КТ. Результаты сравнивались с компенсацией аберраций, рассчитанной по неоднородной модели КТ. Было получено, что максимум давления и эффективность при использовании моделей оКТ и МРТ в этом случае значительно уступают компенсации аберраций по нКТ, так амплитуда в максимуме после коррекции падает до 87 и 65% от случая нКТ, соответственно. Таким образом, модель нКТ с учетом всех неоднородностей тканей предпочтительна для проведения компенсации аберрации в случае реального облучения.

Было проведено также исследование влияния неоднородностей различных тканей на возможности коррекции аберраций. Для этого фазы на излучателе рассчитывались по моделям, в каждой из которых одна из тканей (кожа, череп, мозг) считалась однородной. Было получено, что наибольшее влияние при этом вносят неоднородности черепа. Замена неоднородного сегмента черепа однородным со средними значениями скорости звука и плотности, рассчитанными по конусу облучения, не позволяет точно восстановить положение фокуса после компенсации аберраций и дает наибольшее смещение фокуса относительно геометрического центра кривизны излучателя (2 мм вдоль оси излучателя). Компенсация по такой модели дает также наиболее низкие показатели максимума давления в фокусе (p_A/p₀ = 37.11) и эффективности (13.6%) по сравнению с компенсацией аберраций, рассчитанной по моделям с однородной кожей или мозгом.

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

Работа выполнена при поддержке Программы развития МГУ (проект № 23- Ш06-02) и стипендии фонда развития теоретической физики и математики “Базис” № 22-2-10-6-1.

Авторлар туралы

D. Chupova

Lomonosov Moscow State University

Хат алмасуға жауапты Автор.
Email: daria.chupova@yandex.ru
Ресей, Moscow

P. Rosnitskiy

Lomonosov Moscow State University

Email: daria.chupova@yandex.ru
Ресей, Moscow

O. Solontsov

Lomonosov Moscow State University

Email: daria.chupova@yandex.ru
Ресей, Moscow

L. Gavrilov

Lomonosov Moscow State University

Email: daria.chupova@yandex.ru
Ресей, Moscow

V. Sinitsyn

Lomonosov Moscow State University

Email: daria.chupova@yandex.ru
Ресей, Moscow

E. Mershina

Lomonosov Moscow State University

Email: daria.chupova@yandex.ru
Ресей, Moscow

O. Sapozhnikov

Lomonosov Moscow State University

Email: daria.chupova@yandex.ru
Ресей, Moscow

V. Khokhlova

Lomonosov Moscow State University

Email: daria.chupova@yandex.ru
Ресей, Moscow

Әдебиет тізімі

Қосымша файлдар