Patient-Specific Biomechanical Modeling of the Lung Tumor for Radiation Therapy

Biomechanical patient-specific computational modelling and simulation of the respiratory system is one of the important areas of research in medical imaging and radiation therapy. Lung tumor respiratory motion during irradiation reduces the target coverage and increases dose deposition within healthy tissues. The respiratory motion modifies both the shape and the position of internal organs. This is liable to degrade the quality of radiation treatment of cancer. Breathing is an active and complex process where the respiratory motion is non-reproducible, and the breathing periodicity, amplitude and motion path of patients’ organs, Intra-cycle and inter-cycle variations of respiratory patterns are observed during the respiration (Ehrhardt et al. 2013). This uncertainty on the prgans position makes necessary the development of a strategy for their position prediction. In this paper, we have developed 4D dynamic and nonlinear biomechanical patient-specific computational model of the respiratory system, based on an automatic algorithm of lung pressure and diaphragm force estimation, during a whole respiratory cycle. Using this model, the deformation of the internal organs, are determined by ribs kinematics together with diaphragm motion, correlated with respiratory surrogate signal(s), such as spirometry and/or the displacement of the thorax surface.


Introduction
Biomechanical patient-specific computational modelling and simulation of the respiratory system is one of the important areas of research in medical imaging and radiation therapy. Lung tumor respiratory motion during irradiation reduces the target coverage and increases dose deposition within healthy tissues. The respiratory motion modifies both the shape and the position of internal organs. This is liable to degrade the quality of radiation treatment of cancer. Breathing is an active and complex process where the respiratory motion is non-reproducible, and the breathing periodicity, amplitude and motion path of patients' organs, Intra-cycle and inter-cycle variations of respiratory patterns are observed during the respiration (Ehrhardt et al. 2013). This uncertainty on the prgans position makes necessary the development of a strategy for their position prediction. In this paper, we have developed 4D dynamic and nonlinear biomechanical patient-specific computational model of the respiratory system, based on an automatic algorithm of lung pressure and diaphragm force estimation, during a whole respiratory cycle. Using this model, the deformation of the internal organs, are determined by ribs kinematics together with diaphragm motion, correlated with respiratory surrogate signal(s), such as spirometry and/or the displacement of the thorax surface.

3D segmentation and reconstruction
In this paper, we have developped a methodology based on CAD-modeling to generate a biomechanical patient-specific model of the respiratory system, issued from medical images (Fig. 1). The main motivations of the proposed methodology come from the difficulties in applying the finite element (FE) method in biomechanical analysis due to the complexity associated with creating subject specific anatomical models. For this reason, we have adopted the following computational procedure: (1) 3D automatic and semi-automatic segmentation: The thorax, the lungs and the thoracic flesh have been segmented automatically using gray-level thresholds algorithms. The human diaphragm was segmented manually. Then, we have used the different segmentation masks of the lungs, thorax, the inner thoracic region and the diaphragm to extract the mediastinum structure, (2) Mesh generated: a 3D surface mesh is created for each segmented organ, using the marching cubes algorithm. Then, the surface meshes are rebuilt as a solid using a semi-automatic surface creation with NURBS. Finally, a good quality mesh with four-node tetrahedral elements is generated using ABAQUS packages.

Biomechanical model
To identify the type of the nonlinearities (large stress or large strain), a statistical study on the different strain deformation based on 4D CT images are performed We have obtained an average mean deformation less than 10%. These results confirm that the approach of small strains (with the large stress) may be globally maintained. In our simulation the lungs and diaphragm are considered as isotropic, compressible, with hyperelastic Saint Venant-Kirchhoff model (HSVK), the thorax and mediastinum behavior is considered linear elastic (LE) behavior. For an isotropic elastic or hyper-elastic material, the nonlinear (NL) elastic energy noted W NL can be written as: where E is the Green-Lagrange strain tensor, and are Lame's coefficients. The second Piola-Kirchhoff stress tensor is given by: Therefore, dynamic simulation using FE method, the equation of motion of a vertex v of the organ mesh can be written: where M v , are respectively; the mass computed from Hounsfield densities and damping coefficients of each vertex. The is the neighborhood of vertex. The F int v are the

Conclusions
We have developed an accurate biomechanical patient-specific model of the respiratory system for a whole respiratory cycle. The preliminary results are quite realistic, compared to the ground-truth landmarks. We can observe that the proposed model is able to predict correctly the respiratory motion. Currently, we are working on acceleration of our FE model to find interactively the correlation between the internal organs motion and the external respiratory surrogate signals during treatment. The model could be a potential tool to reduce the treatment margins and consequently to reduce collateral damages to neighboring healthy tissues.
internal forces calculated by FE method and the F v ext are the imposed forces calculated by our developed automatic tuning algorithm based on inverse FE. To solve the dynamic system, the implicit finite scheme has been chosen for more stability.

Boundary conditions (BC)
For the diaphragm, we have applied muscle forces along the radial direction, follwing the direction of muscle fibers. Simple homogeneous Dirichlet boundary condition is applied to the lower part of the diaphragm and the Lagrange multiplier method is used for the contact model. A surface-to-surface contact model is applied to the lung-chest cavity, to simulate the sliding of the lungs. The frictionless contact surfaces are used to simulate the pleural fluid behavior. To simulate the rib cage kinematics, automatic rigid registration algorithm has been developed (Ladjal et al. 2015). The amplitude of the lung pressure and diaphragm forces are patient-specific, and determined at different respiratory states by an optimization framework.

Results and discussion
The proposed patient-specific biomechanical model has been implemented and evaluated on two real patients P1 and P6 from DIR-Lab Dataset (Castillo et al. 2009), with small and large breathing amplitudes. First, we have evaluated the motion estimation accuracy by comparing the simulation results with ground truth (CT images) on 75 landmarks, where the mechanical properties are given in (Fuerst et al. 2015, Ladjal et al.2015. In our simulation, we have obtained an average mean error and standard deviation ( ± SD ) for all ground-truth landmarks: 2.0 ± 1.4 (mm ; Table 1). Moreover, we have compared the lung tumor trajectories identified in 4D CT scan images, with the target trajectories estimated by FE simulation, during the whole breathing cycle, by comparing and calculating the average Hausdorff distance d H between the 3D mesh surface of the segmented tumor (X) and predicted FE lung tumor (Y): The Figure 2 demonstrates that the lung tumor position estimation is accurate with mean error less than 3 mm.