predictive modeling of AAA growth
F. Graph based total variation on non-uniform meshes
The definition of the graph based version of the TV functional
T Vw(b) =X
i
X
j
wG(i, j)(bj−bi)2+2
1/2
, (F.1)
cf. (3.88), allows for a very flexible representation of functions with bounded variation on complex domains. The spatial ‘structure’ of the domain has to be encoded in the adjacency matrix[WG]ij =wG(i, j). Thereby, the term ‘structure’ refers to the spatial connectivity as well as the spatial scale of this connectivity. The former is represented by the connectivity pattern of the adjacency matrix WG, the latter by the values of its entries, the so called edge-weight. For signals being defined on a mesh, this structure is naturally given by the mesh-connectivity, see figure F.1. A representation of signals in
bj
bj bj
bj
bi
bj
bj bj
bj
bi
bj
bj
bj
Figure F.1.:Different connectivity setups for an element-wise parameter layout.
terms of the element-wise basis (2.105) lends itself towards the definition of the vertices of the graph underlying the adjacency matrixWG to be the elements of the mesh. The most straightforward definition of the connectivity among these vertices, i.e., the edges of the graph, is then given by the faces dividing one element from another. At the same time, a face’s area can directly be used as a measure of the scale of connectivity in the definition of the edge-weightswG, see figure F.1 (left). Another possibility of connectivity is given by the definition of edges given by any lower dimensional entity (faces, lines, nodes). However, the natural definition of the spatial scale via a unique entity defining edges is lost. Instead, the edge-weight can be defined via the spatial distance between element centers in this case, see figure F.1 (right). Throughout the work presented in this thesis, the connectivity of the adjacency matrix is defined solely via element faces
since this choice generates graphs that are sparser compared to the graphs defined from connectivity given by faces, lines and nodes.
By denoting the area aij as the area of the face dividing element Ei from element Ej, the edge-weight wG(i, j) can be defined by
wG(i, j) = ( aˆ
aij ⇐⇒ aij 6= 0
0 otherwise, (F.2)
whereby the average areaˆa=Pnf
i=1aiis computed as the average of all nf internal faces ai. This definition is motivated by the regularizing effect of the norm (3.88) on a function defined as element-wise constant function on a mesh. For the convenience of easy meshing or in the case of mesh size adaptivity, meshes are not necessarily homogeneous. But this inhomogeneity should not affect the regularizing properties of the TV norm. In fact, in the sense of the inverse problem set up in chapter 3, the resulting solution should be penalized in a spatially uniform manner by the TV regularization. This is achieved by (F.2), which is shown in the following using synthetic data.
To this end, the synthetic model introduced in chapter 6 is used, see figure 6.1. This model incorporates an inhomogeneous distribution of the mesh size, see figure F.2, due to the structured meshing with hexahedral elements.
Figure F.2.:Local mesh size in terms of element volume.
The synthetic inverse problem is defined in terms of measurements given as surface data and the likelihood given by (6.8). This measured surface is constructed by the solution of the forward problemD=A(θ=cϑ,P)and its respective push-forward action on the outer model surfaceϕA(cϑ,P)(ΓO). In contrast to the problem setup in chapter 6, synthetic noise is applied on this surface in the sense of centered Gaussian noise in the space of currents (σN2 = 0.001 mm2, σW2 = 1.0 mm2), see chapter 3.2. The spatial scale of this noise is not defined by the discretization of the surface but by the spatial scale σW of the kernel defining the space W, see figure 3.7. Despite the change in the size of the triangulation of the model surfaces, the action of the synthetic noise is propagated in a spatially uniform manner via the application of the surface current perspective. Thus, beside the jumps at the patch boundaries, the inverse solution is expected to show a
constant spatial variability, given that the prior penalizes uniformly. This property is verified for the MAP solution of the inverse problem defined by the likelihood in terms of currents and the TV prior (3.93) (αtv = 1.0, = 1.0e−2). For its definition, the use of (F.2) is contrasted with a uniform weight function given by
wG(i, j) =
(1 ⇐⇒ aij 6= 0
0 otherwise. (F.3)
The corresponding MAP solutions are labeledcϑ,Ufor the uniform edge-weights according to (F.3) andcϑ,Afor the area weighted edge-weights according to (F.2), see figure F.3. It
Figure F.3.:MAP solutions cϑ,U corresponding to the uniform edge-weights (left) and cϑ,A
corresponding to the area weighted edge-weights (right). Both solutions reflect the patch-wise character from the reference solution.
can be seen that both solutions reflect the patch-wise character of the ground truthcϑ,P. Visually, a difference in the spatial variability within each of the patches is difficult to detect due to the small local variability compared to the overall range of values covered by the solution. However, this overall range cancels out by analyzing the spatial gradient1
|∇cϑ|, see figure F.4. It can be seen that the spatial gradient, except for the transition between the patches, shows the uniform distribution for the solution obtained with area weighted edge-weights. The solution obtained with uniform edge-weights shows a much higher gradient in regions with a significantly reduced element size.
This example demonstrated that the use of the area weighted edge-weight (F.2) avoids the effects of a heterogeneous mesh size to be reflected in the modeling of prior infor-mation. In contrast, the uniform weight function (F.3) does not account for the spatial variability of smoothness implied by a variation in mesh size. Therefore, definition (F.2) is applied throughout this thesis.
1The gradient is approximated by projecting the element-wise constant cell data to the point data and followingly using the elemental shape functions to interpolate the gradient at the cell center. This gradient is computed by Paraview’s GradientOfUnstructuredDataSet function using vtkCellderiva-tives.
Figure F.4.:Spatial gradient norm|∇cϑ,U|corresponding to the uniform edge-weight (left) and spatial gradient norm |∇cϑ,A| corresponding to the area weighted edge-weights (right). The gradient at the transition between the patches is not shown due to the large magnitude which is outside the chosen color range.
Abbreviations
AAA Abdominal aortic aneurysm BFGS Broyden-Fletcher-Goldfarb-Shanno CI Credible interval
CT Computed tomography
ESS Effective sample size FD Finite difference
FE Finite element
FEM Finite element method FOM Full order model
FSI Fluid structure interaction G&R Growth and remodeling HPD Highest posterior density ILT Intraluminal thrombus LBFGS Limited memory BFGS
LDDMM Large deformation diffeomorphic metric mapping MAP Maximum a posteriori
MC Monte Carlo
MCMC Markov chain Monte Carlo MPI Message passing interface MRF Markov random field
MRI Magnetic resonance imaging
MULF Modified updated Lagrangian formulation OpenMP Open Multi-Processing
PCA Principle component analysis PDE Partial differential equation
PM Posterior mean
POD Proper orthogonal decomposition PTC Pseudotransient continuation RKHS Reproducing kernel Hilbert space RPI Rupture potential index
SEF Strain energy function SMC Sequential Monte Carlo SNR Signal-to-noise ratio
SPD Symmetric positive definite SVD Singular value decomposition
TV Total Variation
UQ Uncertainty quantification US Ultrasound
VB Variational Bayes
Nomenclature
General
(·)> Transpose of(·)
δA[δB] Variation ofA under admissible variationδB
∆A[δB] Variation ofA under admissible variationδB, used ifA=δC