# Robust statistical methods for vibration-based system identification and damage diagnosis

Michael Döhler

<figure><img src="image_1.png" /></figure>

## To cite this version:

Michael Döhler. Robust statistical methods for vibration-based system identification and damage diagnosis. Dynamique, vibrations. Université de Rennes, 2023. tel-04396972

HAL Id: tel-04396972

https://inria.hal.science/tel-04396972

Submitted on 16 Jan 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
HABILITATION À DIRIGER LES RECHERCHES

UNIVERSITÉ DE RENNES

Ecole doctorale Matisse

Domaine : Signal, Image, Vision

présentée par

Michael Döhler

préparée au centre de recherche Inria
de l'Université de Rennes

Robust statistical methods for
vibration-based system identification
and damage diagnosis

Soutenue à Rennes le 7 décembre 2023

devant le jury suivant :

Eleni Chatzi

Professor, ETH Zurich, Switzerland / rapporteure

Spilios Fassois

Professor, University of Patras, Greece / rapporteur

Joseph Morlier

Professor, ISAE-SUPAERO, France / rapporteur

Geert Lombaert

Professor, KU Leuven, Belgium / examinateur

Dmitri Tcherniak

Senior Research Engineer, Hottinger Brüel & Kjær, Denmark / examinateur

Laurent Mevel

Directeur de Recherche, Inria, France / président
<figure><img src="image_3.png" /></figure>
Acknowledgments

First of all, I would like to thank all members of the jury for assessing my work. I am par-
ticularly grateful that Eleni Chatzi (professor at ETH Zurich), Spilios Fassois (professor at
University of Patras) and Joseph Morlier (professor at SUPAERO) have accepted the labo-
rious task of rapporteur of this manuscript. Furthermore, I sincerely thank Geert Lombaert
(professor at KU Leuven) and Dmitri Tcherniak (senior research engineer at HBK) for their
participation as examinateur, and Laurent Mevel for accepting the role as president of the
jury. Thanks for attending the defense in person in Rennes or through video-conference, the
pertinent feedback and inspiring discussions are highly appreciated.

I am grateful to all the academic and industrial collaborators that I had the pleasure
of working with over the years. They have inspired and shaped my research. My deepest
gratitude goes to Laurent Mevel, who advised my as a PhD student and gave me the taste of
research on mathematical methods for structural health monitoring. I am grateful to Palle
Andersen from Structural Vibration Solutions A/S for hosting me half a year during my
PhD studies, and to Dionisio Bernal at Northeastern University and Falk Hille at BAM for
hosting me each for one year as a young postdoc, giving me indispensable background in civil
engineering and applications.

My gratitude goes to Luciano Marin, Delwar Bhuyan, Saeid Allahdadian, Szymon Greś, Alexander Mendler and Eva Viefhues, who will recognize aspects of their PhD work in this manuscript. It has been a pleasure to advise and to work with them. The continued collaboration with Szymon and Alexander after the start of their own academic careers is particularly rewarding.

Thanks also, and especially, to all the members of the I4S research team at Inria in Rennes
and in Nantes, researchers, postdocs and engineers, PhD students, assistants, as well as the
support staff at Inria and Université Gustave Eiffel. It is a pleasure to work together, and
very rewarding to be in this inspiring environment between theory and applications. Thanks
in particular to Laurent Mevel for the great support over the years, his advice, trust in me,
and scientific example.

Last but not least, I wish to thank my family for their never-ending support. Thanks
Mélissande, and thanks Marius and Mathilde, for the love and joy you give me every day.
<figure><img src="image_4.png" /></figure>
# Contents

## Introduction 5

### 1 Subspace-based system identification 9

1.1 Introduction 9

1.2 Efficient multi-order identification 11

1.3 Efficient multi-setup identification 13

1.4 Identification under unknown periodic inputs 16

1.5 Input matrix estimation in combined deterministic-stochastic identification 18

1.6 Conclusions and outlook 20

### 2 Uncertainty quantification in system identification 21

2.1 Introduction 21

2.2 Efficient multi-order computation 23

2.3 Uncertainty in combined deterministic-stochastic identification 25

2.3.1 Matrices A and C 26

2.3.2 Matrices B and D 27

2.4 Uncertainty of modal indicators: second-order framework 28

2.5 Conclusions and outlook 30

### 3 Uncertainty quantification for deterministic damage diagnosis 31

3.1 Introduction 31

3.2 Statistical damage localization with Damage Locating Vector approach 33

3.3 Uncertainties in finite element model updating 36

3.4 Conclusions and outlook 38

### 4 Sensitivity-based statistical damage diagnosis 39

4.1 Introduction 39

4.2 Robustness for local approach methods 41

4.2.1 Introduction 41

4.2.2 Robustness to change in ambient excitation 43

4.2.3 Robustness to uncertainties in reference 45

4.2.4 Robustness to external nuisance 47

4.2.5 Localization: robustness to over-parametrization 50
<table><tr><td>4.2.6 Numerical robustness of test computation</td><td>52</td></tr><tr><td>4.3 Kalman filter-based additive change detection and isolation</td><td>53</td></tr><tr><td>4.4 Conclusions and outlook</td><td>55</td></tr><tr><td><b>5 Damage diagnosis performance and reliability</b></td><td><b>57</b></td></tr><tr><td>5.1 Introduction</td><td>57</td></tr><tr><td>5.2 Detectability and localizability</td><td>58</td></tr><tr><td>5.3 Structural reliability updating and value of information</td><td>61</td></tr><tr><td>5.4 Conclusions and outlook</td><td>64</td></tr><tr><td><b>Perspectives</b></td><td><b>65</b></td></tr><tr><td><b>Publications</b></td><td><b>69</b></td></tr><tr><td><b>References</b></td><td><b>81</b></td></tr></table>
# Introduction

This manuscript summarizes my research activities between 2008 – when I started with my PhD – and 2023 for the French habilitation à diriger des recherches (HDR).

My research is motivated by vibration analysis of mainly civil structures, where I try to bring together the statistical aspects of system identification and fault diagnosis with engineering specifics related to modal analysis and damage diagnosis. I put a strong emphasis on the consideration of statistical uncertainties as well as numerical analysis, with the objectives of developing methods for Operational Modal Analysis (OMA) and Structural Health Monitoring (SHM) that are fast and robust under typical application conditions. In particular, they should ensure the numerical feasibility for realistic problem sizes as well as robustness to noise for structures in operation, such as long-span bridges or high-rise buildings. The goal are automated and robust methods for SHM-based structural diagnosis, on the way to intelligent infrastructure. In my work, theoretical developments are – whenever possible – motivated and illustrated by applications to structures in operation, showing also their effectiveness for industrial problems and leading to various transfers.

The basic problem consists in the analysis of a linear system from measured vibration data of structures in operation, and its monitoring over time, possibly decades, to ensure the serviceability of a structure and to prevent catastrophic failures due to structural damages caused e.g. by corrosion or fatigue. The topic is multidisciplinary between mathematical modelling, statistical signal processing, automatic control and civil engineering, and calls for methods for system identification (to obtain the structural vibration characteristics from data) and fault diagnosis (to detect, locate and quantify changes in the system due to damage from data, in connection with physical models).

Challenges for system identification or fault detection approaches that have originated in the automatic control community are, amongst others, the atypical nature of the considered linear systems: inputs are usually unknown (unmeasured, uncontrolled and non-stationary ambient excitation), the systems are of high model order and the number of outputs (sensors) can be very large. Challenges for traditional modal analysis and monitoring methods that have originated in the civil engineering community are, amongst others, the consideration of statistical uncertainties due to data-based estimates. In this context, my contributions to the field are mainly related to the robustness of methods under realistic application conditions, in particular by considering and evaluating statistical uncertainty, improving numerical efficiency or removing nuisance, with the overall goal of developing automated methods for monitoring with a solid mathematical foundation. My work can be divided into the following
<figure><img src="image_5.png" /><figcaption>Figure 1 – Overview of the chapters and their interaction.</figcaption></figure>

main parts:

**Subspace-based system identification (Chapter 1):** Subspace methods are well-established methods for system identification. They have a strong theoretical background from automatic control, and have very favorable properties for applications related to modal analysis, like being parametric and direct methods, consistency even under non-stationary excitation or asymptotic normality. Critical issues in the context of modal analysis are a large numbers of sensors, unknown model order and noise due to unmeasured ambient excitation. I have contributed with subspace methods that render modal parameter identification computationally more efficient, and adapted subspace methods for typical vibration related engineering problems like for cases with multiple sensor setups or when nuisance from unknown periodic excitation is present.

**Uncertainty quantification in system identification (Chapter 2):** Noise is an issue in all applications dealing with identification from vibration data, so that uncertainty quantification is crucial for evaluating the precision of any result. However, the subspace methods only provide point estimates, but not their uncertainty bounds. Instead, the uncertainty bounds can be obtained with a first-order scheme for sensitivity-based covariance propagation, which is, however, computationally taxing and numerically infeasible on realistic problem sizes in a direct implementation. I have contributed with a mathematical reformulation yielding a fast and memory-efficient computation, which has made the application to vibration data of realistic structures possible. The developed methodological framework has opened the field for uncertainty quantification of other methods, where I have contributed with uncertainty quantification for a range of output-only and input/output subspace identification methods, as well as for the classical modal indicators MAC (Modal Assurance Criterion) and MPC (Modal Phase Collinearity).

**Uncertainty quantification for damage diagnosis (Chapter 3):** System identification from vibration measurements is a fundamental task for many engineering applications. Amongst others, the obtained modal parameter estimates (or related quantities) are used
in damage localization or quantification methods. Such methods are usually based on me-
chanical or physical grounds and use the modal parameters in a deterministic way, but they
rarely consider data-based uncertainties when inferring about the damage location or extent.
After the developments in the previous chapter, we can profit from the available modal pa-
rameter uncertainties, with the goal of enhancing the damage diagnosis methods. I have
contributed with the integration of estimation uncertainties in the well-established damage
locating vector (DLV) framework, in order to develop a robust statistical evaluation for dam-
age localization with this method class. I also have contributed with the consideration of
data-based uncertainties in Finite Element (FE) model updating, by evaluating and integrat-
ing the uncertainties in the optimization procedure to conceive an efficient updating method
for damage quantification.

**Sensitivity-based statistical damage diagnosis (Chapter 4):** In this chapter, the damage diagnosis problem is considered in a statistical framework from the outset, where a damage-sensitive feature is computed from measurement data and analyzed by statistical hypothesis testing. The link to physical system properties is made by the feature's sensitivities with respect to a physical parameterization. This allows noise robustness in the feature evaluation on the one hand, and a simple link between data and physical parameters for damage diagnosis without the need of model updating on the other hand. In contrast, the consideration of statistical uncertainties and physical system properties was made the other way round in the previous chapter, where data-related statistical uncertainties were included in damage diagnosis methods that were already based on physical properties. I contributed to damage detection and localization within the well-established *asymptotic local approach framework* by developing theory for typical application cases, ensuring robustness for real applications. Developments include robustness towards naturally changing excitation properties, consideration of uncertainties in the data-driven reference model and of nuisance parameters like temperature changes, as well as robustness for damage localization via sensitivity-based clustering approaches. Finally, we have developed a new sensitivity-based Kalman filter-based framework for statistical damage diagnosis as an alternative to the previous local approach framework.

**Damage diagnosis performance and reliability (Chapter 5):** An important feature of the statistical frameworks for damage diagnosis in the previous chapter is the characterization of the probability distributions of the damage detection and localization tests for both reference and damaged states. On this basis, we can define and analyze the *damage detectability and localizability* of the respective methods under the typical data-based uncertainties. This way, the performance of the considered damage diagnosis method can be assessed even before damage occurs, which is an essential information for stakeholders or decision makers to deploy SHM in practice. This information also allows to establish a link between SHM-based damage diagnosis and structural reliability and performance assessment in connection with structural deterioration models in a Bayesian setting. I have contributed to data-based updating of the structural reliability and SHM system performance analysis, as well as an analysis of the value of monitoring information.
# Subspace-based system identification

## 1.1 Introduction

The identification of dynamic system characteristics from vibration measurements is a fundamental task in engineering. Amongst others, subspace-based system identification methods are well-suited for this purpose, identifying the system matrices of a linear time-invariant (LTI) state-space model that describes the dynamic system behavior [1]. Indeed, when the vibration behavior of a structure can be assumed as linear, it can be modelled by the equation of motion

$$ \mathcal{M}\ddot{z}(t) + \mathcal{C}\dot{z}(t) + \mathcal{K}z(t) = f(t), \quad (1.1) $$

where $\mathcal{M}, \mathcal{C}$, and $\mathcal{K} \in \mathbb{R}^{m \times m}$ are the mass, damping, and stiffness matrices, with $m$ denoting the degrees of freedom and $t$ the continuous time. Vector $z \in \mathbb{R}^m$ contains the displacements at the degrees of freedom of the structure, and $f(t) \in \mathbb{R}^m$ is the excitation force. When sampling system (1.1) at a discrete time step $\tau$, measuring the system outputs with $r$ sensors and transforming it to a first-order system, a discrete-time state-space representation [2] is conveniently obtained, with the advantage that the well-known system identification methods originating from the automatic control community can be applied in this format.

In the context of *operational modal analysis* (OMA), where the monitored structures are excited by natural ambient disturbances (road traffic, wind, sea water wave, etc.) that are unknown in most practical situations, $f(t)$ is unknown and the corresponding *output-only* LTI state-space model writes

$$ \begin{cases} x_{k+1} = Ax_k + w_k \\ y_k = Cx_k + v_k \end{cases} \quad (1.2) $$
where $x_k \in \mathbb{R}^n$ is the internal state vector, $y_k \in \mathbb{R}^r$ is the output vector corresponding to available sensor measurements (displacement, velocity, and/or acceleration), $w_k \in \mathbb{R}^n$ and $v_k \in \mathbb{R}^r$ are random noises modeling natural excitation of the structure and measurement noise, $A \in \mathbb{R}^{n \times n}$ is the state transition matrix and $C \in \mathbb{R}^{r \times n}$ is the output matrix. Besides the unknown ambient excitation forces acting on a structure, some of the input forces are available in different engineering applications, such as shaker tests of bridges [3], structural health monitoring of wind turbine blades [4], aircraft control [5] and various laboratory testing applications. In this context of operational modal analysis with exogenous inputs (OMAX), some parts of $f(t)$ are known and the corresponding input/output LTI state-space model writes

$$ \left\{ \begin{aligned} x_{k+1} &= Ax_k + Bu_k + w_k \\ y_k &= Cx_k + Du_k + v_k \end{aligned} \right. \qquad (1.3) $$

where $u_k \in \mathbb{R}^u$ is the known input, $B \in \mathbb{R}^{n \times u}$ and $D \in \mathbb{R}^{r \times u}$ are the input and feedthrough matrices.

It is our goal to estimate matrices $(A, C)$ from measurement data to obtain the modal parameters (natural frequencies, damping ratios, mode shapes) based on the eigenvalues and eigenvectors of $A$ and $C$. If inputs are available, matrices $(B, D)$ can be estimated in addition to infer about the complete input-output behavior of the structure, e.g., by its transfer function. For these tasks, several stochastic and combined deterministic-stochastic subspace methods are available [1], [6]–[9].

While these methods originate from the automatic control community, the problem of vibration analysis is quite untypical therein, since the model orders $n$ are unknown and can be large, the number of sensors can be large (or sensors are moved between different measurement setups) and the identification can be perturbed by non-white excitation noises or periodic excitation leading to additional spurious noise modes or harmonic modes in addition to the structural modes. In this chapter, we have considered the problems of subspace-based system identification at high model orders with many sensors or multiple sensor setups, identification under periodic excitation and identification with (partly) known inputs.

## Contributions

My contributions are mainly related to numerical efficiency and robustness of subspace methods under realistic application conditions:

* Numerically efficient identification at multiple model orders. The method significantly speeds up the computation of the classic stabilization diagram for modal analysis, in particular for high model orders (Section 1.2).

*Context: This work was part of my PhD, supervised by Laurent Mevel, at Inria.*

* Fast and scalable subspace identification using multiple sensor setups. The method identifies a global model based on the different setups, and takes changing excitation properties between setups and unknown model order into account (Section 1.3).

*Context: This work was also part of my PhD, supervised by Laurent Mevel.*
* Robustness of subspace identification under unknown periodic inputs. Such inputs stem e.g. from periodic movement of rotating machinery, which perturb modes close to the rotation frequency, and the method rejects their effect on the measurements (Section 1.4).

Context: Collaboration with Szymon Greś during his PhD at Aalborg University [10], and with Laurent Mevel.

* Consistency of estimation of $(B, D)$ for a full input-output characterization. For the case when input measurements are actually available, we have proven consistency of $(B, D)$ in a more theoretical exercise (Section 1.5), before analyzing their uncertainty in the subsequent chapter.

Context: Postdoc of Szymon Greś at Inria, together with Laurent Mevel.

## 1.2 Efficient multi-order identification

The stabilization diagram is a standard tool in operational modal analysis, where the model order is unknown and the noise terms related to the unknown ambient excitation are in general not white in models (1.2) or (1.3). Under the assumption that (spurious) modes related to the noise vary when identified at different model orders, while the physical structural modes remain stable, the stabilization diagram allows to distinguish both [7]. This requires their identification at multiple, possibly over-specified model orders, which can be an important computational burden and a problem for real-time monitoring or in environments with limited computational capacities like in embedded sensor platforms. Previously, literature on computational efficiency for subspace methods concerned the reduction of the processed data by using projection channels [7], fast processing of the measurement data [11], [12] or recursive algorithms [13], but not the identification at multiple model orders.

In subspace identification, matrix $A$ is obtained as the solution of a least-squares problem from the observability matrix estimate, which is obtained from a singular value decomposition of a matrix $\hat{\mathcal{H}}$ and its truncation at the model order

$$ \hat{\mathcal{H}} = [\hat{U}_1 \quad \hat{U}_2] \begin{bmatrix} \hat{\Sigma}_1 & 0 \\ 0 & \hat{\Sigma}_2 \end{bmatrix} \begin{bmatrix} \hat{V}_1^T \\ \hat{V}_2^T \end{bmatrix}, \quad \hat{\mathcal{O}} = \hat{U}_1 \hat{\Sigma}_1^{1/2}. \qquad (1.4) $$

The construction of matrix $\hat{\mathcal{H}}$ from output data (as well as input data, if available) data Hankel matrices such that the observability matrix defines a basis of its column space is common to all subspace methods and depends on the chosen method, see e.g. details in [1], [7], [8]. Choosing different model orders corresponds to different truncations of the SVD in (1.4), and thus to appending columns to $\hat{\mathcal{O}}$. Since $\hat{A}$ is computed from $\hat{\mathcal{O}}$ at different orders, the key to reducing the computational complexity is the exploitation of the structure of the corresponding least-squares problems.

The observability matrix and its versions without the last and first block rows are respec-
tively defined as

$$ \mathcal{O} = \begin{bmatrix} C \\ CA \\ \vdots \\ CA^p \end{bmatrix}, \quad \mathcal{O}^{\uparrow} = \begin{bmatrix} C \\ CA \\ \vdots \\ CA^{p-1} \end{bmatrix}, \quad \mathcal{O}^{\downarrow} = \begin{bmatrix} CA \\ \vdots \\ CA^p \end{bmatrix}, \qquad (1.5) $$

and the least squares problem for A writes

$$ \mathcal{O}^{\uparrow} A = \mathcal{O}^{\downarrow}. \qquad (1.6) $$

Noting that the observability matrix is supposed to have full column rank, the solution is usually given by $A = (\mathcal{O}^{\uparrow})^{\dagger}\mathcal{O}^{\downarrow}$. This solution needs to be computed for all the model orders $n = 1, \dots, n_{\max}$ in the stabilization diagram. To denote the model order, the index $n$ is added to the matrices $O_n, A_n$ and $C_n$. The main idea to reduce the computational complexity is to perform the computations at the maximal model order $n_{\max}$ first, where the QR decomposition is used to solve the least squares problem instead of the pseudoinverse. With this, the main computational burden is already done, since we have shown in [J2] that the solution of the least squares problem for $n < n_{\max}$ follows nearly effortlessly based on the matrices computed at $n_{\max}$:

**Theorem 1.1** Let the thin QR decomposition of $O_{n_{\max}}^{\uparrow}$ be given with

$$ O_{n_{\max}}^{\uparrow} = Q_{n_{\max}} R_{n_{\max}}, \quad S_{n_{\max}} = Q_{n_{\max}}^{T} O_{n_{\max}}^{\downarrow}, \qquad (1.7) $$

such that $A_{n_{\max}} = R_{n_{\max}}^{-1} S_{n_{\max}}$ is the least squares solution of (1.6) at order $n_{\max}$. Let $n < n_{\max}$, and let $R_{n_{\max}}$ and $S_{n_{\max}}$ be partitioned into

$$ R_{n_{\max}} = \begin{bmatrix} R_n^{(11)} & R_n^{(12)} \\ 0 & R_n^{(22)} \end{bmatrix}, \quad S_{n_{\max}} = \begin{bmatrix} S_n^{(11)} & S_n^{(12)} \\ S_n^{(21)} & S_n^{(22)} \end{bmatrix}, \qquad (1.8) $$

where $R_n^{(11)}, S_n^{(11)} \in \mathbb{R}^{n \times n}$. Then, the state transition matrix $A_n$ at model order $n$, which is the least squares solution of (1.6) at order $n$, satisfies

$$ A_n = (R_n^{(11)})^{-1} S_n^{(11)}. \qquad (1.9) $$

Hence, $A_n$ can be simply obtained based on submatrices of the matrices that are already computed at $n_{\max}$. Thanks to this property, the burden for computing the system matrices for a whole stabilization diagram is essentially reduced to computing the system matrices at only one model order $n_{\max}$. Computing (1.9) for $n = 1, \dots, n_{\max}$ is not very expensive in practice, but still takes $O(n_{\max}^4)$ flops, which can be reduced to $O(n_{\max}^3)$ flops when computing $A_n$ iteratively based on (1.8)–(1.9). The computation times for the conventional computation of the system matrices (solving (1.6) for each model order), the fast computation (based on Theorem 1.1) and the iterative fast computation (computing $A_n$ in (1.9) iteratively) are depicted in Fig. 1.1, showing that the new algorithms are about 100 times faster for realistic problem sizes.
<figure><img src="image_6.png" /><figcaption>Figure 1.1 – Computation times for multi-order identification of system matrices $A_n$ and $C_n$, $n = 1, \dots, n_{max}$ from $O_{n_{max}}$ at different maximal model orders $n_{max}$.</figcaption></figure>

[J2] M. Döhler and L. Mevel, “Fast multi-order computation of system matrices in subspace-based system identification,” *Control Engineering Practice*, vol. 20, no. 9, 882–894, 2012

[CA3] M. Döhler and L. Mevel, “Fast multi-order stochastic subspace identification,” in *Proc. 18th IFAC World Congress*, Milan, Italy, 2011

[CS20] M. Döhler, P. Andersen, and L. Mevel, “Operational modal analysis using a fast stochastic subspace identification method,” in *Proc. 30th International Modal Analysis Conference (IMAC)*, Jacksonville, FL, USA, 2012

[CS32] V. Le Cam, M. Döhler, M. Le Pen, and L. Mevel, “Embedded modal analysis algorithms on the smart wireless sensor platform PEGASE,” in *Proc. 9th International Workshop on Structural Health Monitoring (IWSHM)*, Stanford, CA, USA, 2013

## 1.3 Efficient multi-setup identification

To obtain vibration measurements at many coordinates of a structure with only a few sensors, it is common practice to use multiple sensor setups for the measurements, where some of the sensors (reference sensors) remain fixed while others (moving sensors) are moved from setup to setup. By fusing in some way the corresponding data, this allows to perform modal identification as if there was a very large number of sensors, even in the range of a few hundreds or thousands. A classic approach is individual modal analysis on the data of each setup, and then merging the results by tracking the modes throughout all setups, averaging frequencies and damping ratios, and gluing the mode shape components of the different setups together after normalization with respect to the reference components [14]. However, the mode track-
ing can be cumbersome when some of the modes are weakly excited in some of the setups, or when modes are closely spaced. We pursue a global approach, where the data is described by a global model (since the same system is observed, but at different output coordinates in the different setups). In [15], [16], a method was proposed for covariance-driven subspace identification to merge and to normalize the data from all the setups, before doing the global system identification on it in one step. In this contribution, we have generalized it to a large range of subspace methods, including data-driven subspace methods [1]. With this merging strategy, a theoretically sound and simple identification of the global modal parameters is possible in one run. Furthermore, the merging approach can be used in an iterative scheme where no global Hankel or observability matrix is necessary, but the desired system matrices can be computed iteratively setup by setup, which is a computational advantage if the number of setups is large.

The method is outlined as follows. Instead of a single output record {$y_k$} of the system (1.2), $N_s$ records

$$ \underbrace{\begin{Bmatrix} y_k^{(1,\text{ref})} \\ y_k^{(1,\text{mov})} \end{Bmatrix}}_{\text{Record 1}} \underbrace{\begin{Bmatrix} y_k^{(2,\text{ref})} \\ y_k^{(2,\text{mov})} \end{Bmatrix}}_{\text{Record 2}} \cdots \underbrace{\begin{Bmatrix} y_k^{(N_s,\text{ref})} \\ y_k^{(N_s,\text{mov})} \end{Bmatrix}}_{\text{Record } N_s} \qquad (1.10) $$

are now available collected successively. Each record $j$ contains data $y_k^{(j,\text{ref})}$ of dimension $r^{(\text{ref})}$ from a fixed reference sensor pool, and data $y_k^{(j,\text{mov})}$ of dimension $r^{(j)}$ from a moving sensor pool. To each record $j = 1, \dots, N_s$ corresponds a state-space model in the form

$$ \left\{ \begin{aligned} x_{k+1}^{(j)} &= A x_k^{(j)} + w_k^{(j)} \\ y_k^{(j,\text{ref})} &= C^{(\text{ref})} x_k^{(j)} + v_k^{(j,\text{ref})} && (\text{reference pool}) \\ y_k^{(j,\text{mov})} &= C^{(j,\text{mov})} x_k^{(j)} + v_k^{(j,\text{mov})} && (\text{sensor pool } j) \end{aligned} \right. \qquad (1.11) $$

with a single state transition matrix $A$, since the same system is being observed. As such, all the models (1.11) for $j = 1, \dots, N_s$ are defined in the same state-space basis, to ensure the same matrix $A$ and compatible observation matrices for all setups. The observation matrix $C^{(\text{ref})}$ with respect to the reference sensors is independent of the measurement setup as the reference sensors are fixed throughout the measurements, while the observation matrices $C^{(j,\text{mov})}$ correspond to the moving sensor pool of each setup $j$.

In order to identify system (1.11) on the basis of the measurements (1.10), the basic subspace algorithm is adapted that is based on the estimation of a matrix $\mathcal{H}$ from data (specific to the selected subspace method), from which the observability matrix $\mathcal{O}$ is obtained from its column space (see (1.4)-(1.5)). While matrix $\mathcal{H}$ has such a factorization property $\mathcal{H} = \mathcal{O}\mathcal{Z}$ (with a matrix $\mathcal{Z}$ specific to the selected subspace method), “local” matrices $\mathcal{H}^{(j)}$ computed on each of the $N_s$ data records in (1.10) also yield such a factorization property with different observability matrices $\mathcal{O}^{(j)}$ on the left (due to the different sensor setups), but in general also with different matrices $\mathcal{Z}^{(j)}$ on the right side that prevents a straightforward merging. In particular, obstacles are:

* In general, the matrices $\mathcal{Z}^{(j)}$ depend on different sensor sets. This can be solved with reference-based subspace identification [7].
* For many subspace methods, the matrices $\mathcal{Z}^{(j)}$ depend directly on the data (and thus on different initial states $x_0^{(j)}$), or they are not uniquely defined due to LQ decompositions.

* Matrix $\mathcal{Z}^{(j)}$ is influenced by the noise properties of the system, which may change from setup to setup since multi-setup measurements take place over a longer time period than a single measurement.

* Estimates of matrices $\mathcal{O}^{(j)}$ and $\mathcal{Z}^{(j)}$ for different setups are in general obtained in different state-space bases, so that in fact $\hat{\mathcal{O}}^{(j)}T_j$ and $T_j^{-1}\hat{\mathcal{Z}}^{(j)}$ are identified with an invertible unknown change of basis matrix $T_j \in \mathbb{R}^{n \times n}$ that is different for each setup.

In short, the challenge is to obtain the parts of the observability matrices from the different setups in the same basis, then their merging is straightforward and allows a global system identification, which is moreover generic for any subspace method. This is possible thanks to the reference sensors. For this, the reference and moving sensor parts $\mathcal{O}^{(j,ref)}$ and $\mathcal{O}^{(j,mov)}$, respectively, are selected from the block rows of $\mathcal{O}^{(j)}$. These matrices are in general obtained in different bases, i.e., an unknown change of basis matrix $T_j$ assures that $\mathcal{O}^{(j)}T_j$, $j = 1, \dots, N_s$, are compatible for different setups. Among the setups, one reference matrix $\mathcal{O}^{(\text{ref})} = \mathcal{O}^{(j,ref)}$ is selected. Since all matrices $\mathcal{O}^{(j,ref)}T_j$ are equal, the transformation

$$ \bar{\mathcal{O}}^{(j,mov)} \stackrel{\text{def}}{=} \mathcal{O}^{(j,mov)} \mathcal{O}^{(j,ref)} \dagger \mathcal{O}^{(\text{ref})} \quad (1.12) $$

leads to observability matrix parts in the same basis as $\mathcal{O}^{(\text{ref})}$ that are compatible for merging and consequently global system identification. The resulting method is modular, since the normalization (1.12) is carried out individually on the different setups. It allows an iterative computation for cases with many measurement setups, including an iterative computation for A and the global observation matrix C with respect to the setups. Therefore, this strategy scales well with the number of measurement setups.

To illustrate an application, the mode shapes of Z24 Bridge are shown in Fig. 1.2, which are identified with the described method from vibration data at 251 coordinates in nine different sensor setups.

<figure><img src="image_7.png" /><figcaption>Figure 1.2 – Mode shapes of Z24 Bridge.</figcaption></figure>
[J3] M. Döhler and L. Mevel, “Modular subspace-based system identification from multi-setup measurements,” *IEEE Transactions on Automatic Control*, vol. 57, no. 11, 2951–2956, 2012

[J4] M. Döhler, X.-B. Lam, and L. Mevel, “Uncertainty quantification for modal parameters from stochastic subspace identification on multi-setup measurements,” *Mechanical Systems and Signal Processing*, vol. 36, no. 2, 562–581, 2013

[CS5] M. Döhler, P. Andersen, and L. Mevel, “Data merging for multi-setup operational modal analysis with data-driven SSI,” in *Proc. 28th International Modal Analysis Conference (IMAC)*, Jacksonville, FL, USA, 2010

[CS6] M. Döhler, E. Reynders, F. Magalhães, L. Mevel, G. De Roeck, and Á. Cunha, “Pre- and post-identification merging for multi-setup OMA with covariance-driven SSI,” in *Proc. 28th International Modal Analysis Conference (IMAC)*, Jacksonville, FL, USA, 2010

[CS7] M. Döhler, X.-B. Lam, and L. Mevel, “Crystal clear data fusion in subspace system identification and damage detection,” in *Proc. 5th International Conference on Bridge Maintenance, Safety and Management (IABMAS)*, Philadelphia, PA, USA, 2010

[CS12] M. Döhler, B. Jaishi, L. Mevel, and J. Brownjohn, “Data fusion for system identification of the Humber Bridge,” in *Proc. 29th International Modal Analysis Conference (IMAC)*, Jacksonville, FL, USA, 2011

## 1.4 Identification under unknown periodic inputs

Periodic excitation, for example due to rotating machinery, may render modal parameter estimation difficult. Modes due to the periodic excitation may disturb the estimation of close structural modes, or mask them due to their high energy content. To facilitate the modal parameter estimation, it is desirable to separate the periodic subsignal from the random response signal, and subsequently the modal parameters are estimated with classical methods.

For example, time-synchronous-averaging (TSA) is a method extracting periodic wave-forms from signals by averaging their blocks synchronized in the angular domain. For OMA, this averaged signal is subtracted from the raw measurements, which results in the removal of the periodic frequencies selected to synchronize the blocks [17]. Angle matching is often achieved with tachometer measurements, which is not practical in real-life applications and was attempted to be overcome in the context of TSA in [18]. A family of methods that does not require tachometer measurements is based on the cepstrum, which is an inverse Fourier transform of the logarithm of spectrum [19], [20]. Although the cepstrum is capable to filter the periodic frequency components out of the spectra of the output data, its empirical premise does not ensure the consistency of the resulting modal parameter estimates.

The aim of this work is to develop a consistent subspace identification method for the
estimation of the structural modes while rejecting the influence of the periodic input. Assume that an unknown deterministic periodic force $u(t)$ acts on the system in addition to the random noise input $w(t)$. Denote the resulting states as $x_k^{sys}$, containing the displacements and velocities at the DOFs of the structure under both unknown noise and periodic forces. Then, the continuous-time state space model (see also Eq. (1.3)) writes

$$ \dot{x}^{\text{sys}}(t) = A_c^{\text{sys}}x^{\text{sys}}(t) + \mathbf{b}u(t) + w(t), \quad (1.13) $$

$$ y(t) = C^{\text{sys}}x^{\text{sys}}(t) + \mathbf{d}u(t) + v(t), \quad (1.14) $$

with appropriate matrices $\mathbf{b}$ and $\mathbf{d}$. The periodic force can be modelled as

$$ u(t) = \sum_{i=1}^{h} a_i \sin(\omega_i t + g_i), \quad (1.15) $$

and can be eliminated in model (1.13)–(1.14) by augmenting the state vector with the input dynamics to

$$ \begin{bmatrix} \dot{x}^{\text{sys}}(t) \\ \dot{x}^{\text{per}}(t) \end{bmatrix} = \begin{bmatrix} A_c^{\text{sys}} & A_c^{\mathbf{b}} \\ 0 & A_c^{\text{per}} \end{bmatrix} \begin{bmatrix} x^{\text{sys}}(t) \\ x^{\text{per}}(t) \end{bmatrix} + \begin{bmatrix} w(t) \\ 0 \end{bmatrix}, \quad (1.16) $$

$$ y(t) = \begin{bmatrix} C^{\text{sys}} & C^{\text{per}} \end{bmatrix} \begin{bmatrix} x^{\text{sys}}(t) \\ x^{\text{per}}(t) \end{bmatrix} + v(t), \quad (1.17) $$

where $A_c^{\text{per}}$ relates to the input dynamics whose discrete-time counterpart has all eigenvalues on the unit circle. Hence, the measurements can be described by a state-space model where the periodic input dynamics appear as a subsystem in addition to the structural system of interest, and subspace identification is still consistent [21]. In particular, it allows the estimation of the periodic modes in a preliminary step, and the corresponding state estimates with a non-steady state Kalman filter, i.e., the part of the states related to $x^{\text{per}}$ (in discrete time). They are then removed from the original output signal by an orthogonal projection that is adapted in the context of subspace methods. Consequently, the data is “cleaned” from

<figure><img src="image_8.png" /><figcaption>Figure 1.3 – Stabilization diagram of the natural frequency estimates from measurements of a plate with periodic excitation before (left) and after (right) the removal of the periodic subsignal.</figcaption></figure>
the periodic contributions, and we have proved that the modes of the structural system are consistently estimated. This allows a better estimation of structural modes close to periodic ones, and an easier interpretation of standard data analysis procedures, like the stabilization diagram. An illustration on the identification from vibration data of a plate with periodic excitation is given in Fig. 1.3.

[J19] S. Greś, M. Döhler, P. Andersen, and L. Mevel, “Kalman filter-based subspace identification for operational modal analysis under unmeasured periodic excitation,” *Mechanical Systems and Signal Processing*, vol. 146, 206996, 2021

## 1.5 Input matrix estimation in combined deterministic-stochastic identification

For the identification of modal parameters, only the matrices (A, C) of the state-space model are required, which can be obtained by output-only measurements under ambient excitation. Besides the unknown ambient excitation forces acting on a structure, some of the input forces are available in different engineering applications, such as shaker tests of bridges [3], structural health monitoring of wind turbine blades [4], aircraft control [5], [22] and various laboratory testing applications. Then, the knowledge of inputs enhances the identification performance when using input/output system identification methods – also called combined deterministic-stochastic identification methods. Moreover, the complete set of system matrices, including matrices (B, D) of system (1.3) can be identified. These matrices are required when evaluating the full dynamic system properties, e.g., for computing the transfer function of the system for damage detection [23], damage localization [24] model updating [25], and modal scaling [9].

In the context of subspace methods, several ways of identifying (B, D) were described in [1], whose implementation can be quite complex. A simpler method in the subspace context was sketched in [26], which is independent from the method that identifies (A, C) and is therefore applicable in conjunction with any subspace identification method. Moreover, uncertainty quantification of this method is conveniently possible (see Section 2.3.2), which was the main motivation of this line of work. The detailed development of the identification method and in particular its proof of consistency (also as a prerequisite for the uncertainty quantification) deserve their own attention in this section.

The identification of (B, D) is based on iterating the state-space model (1.3) for *p* time steps, leading to the relation

$$ \mathcal{Y}^{+} = \Gamma \mathcal{X}^{+} + \overline{H} \mathcal{U}^{+} + \mathcal{V}^{+}, \qquad (1.18) $$

where

$$ \mathcal{Y}^{+} = \frac{1}{\sqrt{N}} \begin{bmatrix} y_1 & y_2 & \cdots & y_N \\ y_2 & y_3 & \cdots & y_{N+1} \\ \vdots & \vdots & \ddots & \vdots \\ y_{p+1} & y_{p+2} & \cdots & y_{N+p} \end{bmatrix}, \quad \mathcal{X}^{+} = \frac{1}{\sqrt{N}} [x_1 \ x_2 \ \cdots \ x_N], \qquad (1.19) $$
$$ \overline{H} = \begin{bmatrix} D & 0 & \dots & 0 \\ CB & D & \dots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ CA^{p-1}B & CA^{p-2}B & \dots & D \end{bmatrix}, \quad U^+ = \frac{1}{\sqrt{N}} \begin{bmatrix} u_1 & u_2 & \dots & u_N \\ u_2 & u_3 & \dots & u_{N+1} \\ \vdots & \vdots & \ddots & \vdots \\ u_{p+1} & u_{p+2} & \dots & u_{N+p} \end{bmatrix}, \quad (1.20) $$

and $V^+$ collects the noise terms. Multiplying (1.18) with left null space estimate $\hat{U}_2^T$ (like in (1.4)) from the left and with $U^+^T$ from the right yields

$$ \hat{U}_2^T \hat{R}_1 = \underbrace{\hat{U}_2^T \Gamma}_{\to 0} \mathcal{X}^+ \mathcal{U}^{+T} + \hat{U}_2^T \overline{H} \hat{R}_2 + \hat{U}_2^T \underbrace{\mathcal{V}^+ \mathcal{U}^{+T}}_{\to 0}, \quad (1.21) $$

$$ \hat{U}_2^T \hat{R}_1 \hat{R}_2^{-1} = \hat{U}_2^T \overline{H} + o(1) \quad (1.22) $$

where $\hat{R}_1 = \mathcal{Y}^+ \mathcal{U}^{+T}$, $\hat{R}_2 = \mathcal{U}^+ \mathcal{U}^{+T}$, and $o(1)$ is a matrix whose norm converges to zero if $N \to \infty$. Then, matrix $B$ and $D$ in matrix $\overline{H}$ can be estimated by rearranging Eq. (1.22): Partitioning $\hat{\mathcal{M}} = \hat{U}_2^T \hat{R}_1 \hat{R}_2^{-1} = [\hat{\mathcal{M}}_1 \dots \hat{\mathcal{M}}_{p+1}]$ and $\hat{U}_2^T = [\hat{\mathcal{L}}_1 \dots \hat{\mathcal{L}}_{p+1}]$ leads to

$$ \underbrace{\begin{bmatrix} \hat{\mathcal{M}}_1 \\ \hat{\mathcal{M}}_2 \\ \vdots \\ \hat{\mathcal{M}}_{p+1} \end{bmatrix}}_{=\hat{\mathcal{M}}^v} = \underbrace{\begin{bmatrix} \hat{\mathcal{L}}_1 & \hat{\mathcal{L}}_2 & \dots & \hat{\mathcal{L}}_p & \hat{\mathcal{L}}_{p+1} \\ \hat{\mathcal{L}}_2 & \hat{\mathcal{L}}_3 & \dots & \hat{\mathcal{L}}_{p+1} & 0 \\ \vdots & \vdots & \ddots & \vdots & \vdots \\ \hat{\mathcal{L}}_{p+1} & 0 & \dots & 0 & 0 \end{bmatrix}}_{=\hat{\mathcal{L}}} \underbrace{\begin{bmatrix} I_r & 0 \\ 0 & \Gamma \end{bmatrix}}_{=O_s} \begin{bmatrix} D \\ B \end{bmatrix} + o(1), \quad (1.23) $$

from where least-squares estimates of $(B, D)$ are obtained by

$$ \begin{bmatrix} \hat{D} \\ \hat{B} \end{bmatrix} = (\hat{L}\hat{O}_s)^{\dagger} \hat{\mathcal{M}}^v. \quad (1.24) $$

To prove consistency, the consistency of $\hat{R}_1$, $\hat{R}_2$ and $\hat{\Gamma}$ is required, which is given under standard assumptions (like inputs are quasi-stationary, see also [8] for details). Matrix $\hat{U}_2$ is derived from $\hat{\Gamma}$ by the SVD and therefore is a consistent left null space estimate (following from [27]), leading to consistency of $\hat{\mathcal{M}}^v$ and $\hat{L}$. The main contribution of this work is the proof that $L$ has full column rank, which hinges on $\mathcal{L}_{p+1}$ having full column rank. This can be shown by construction of the left null space matrix $U_2$ and follows from the full column rank of $\Gamma$. The details of this proof are found in [J25, Appendix A]. Then, the product $LO_s$ has full column rank, since $O_s$ has full column rank, and the consistency of the estimates of $(B, D)$ in (1.24) follows.

[J25] S. Greś, M. Döhler, N.-J. Jacobsen, and L. Mevel, “Uncertainty quantification of input matrices and transfer function in input/output subspace system identification,” *Mechanical Systems and Signal Processing*, vol. 167, 108581, 2022

[CA14] S. Greś, M. Döhler, and L. Mevel, “Variance computation for system matrices and transfer function from input/output subspace system identification,” in *Proc. 21st IFAC World Congress*, 2020
1.6 Conclusions and outlook

Subspace methods are well-established, having a strong theoretical background from auto-
matic control, and have very favorable properties for applications (like parametric and direct
identification, consistency even under non-stationary excitation, asymptotic normality, ...).
The developments in this chapter render modal parameter identification computationally
more efficient, and adapt subspace methods for typical vibration related engineering prob-
lems, like for cases with multiple sensor setups or when nuisance from unknown periodic
excitation is present. The made developments increase the applicability and robustness of
subspace methods for realistic vibration-based identification problems. In particular, the
efficient multi-order identification (Section 1.2), the multi-setup algorithm (Section 1.3) as
well as the rejection of the periodic excitation (Section 1.4) have been transferred to the
commercial software ARTEMIS Modal [28]. Furthermore, some of these developments have
laid out the basis for efficient uncertainty quantification in Chapter 2.

With the advent of new sensor technologies, such as video camera-based full field dis-
placement or velocity measurements, the number of measured outputs is quickly increasing,
which is a challenge for the subspace methods. Ongoing developments concern the efficient
processing of sensor data in this case. With postdoc Boualem Merainani [J30], [CS69] and in
the ongoing PhD of Zhilei Luo we investigate PCA-based reduction techniques to keep the
computational burden as well as the memory usage feasible, while preserving the full spatial
resolution of the mode shapes [CA18], [CS81].
# Uncertainty quantification in
system identification

## 2.1 Introduction

All estimates from system identification are inherently afflicted with statistical uncertainties due to measurement noise, (at least partly) unknown inputs related to ambient excitation, and limited data length. The quantification of these uncertainties is important for diverse engineering applications related to Operational Modal Analysis or Structural Health Monitoring, e.g., to evaluate confidence bounds of estimates, or to evaluate if changes during monitoring are statistically significant for damage diagnosis.

The uncertainty quantification is inherently tied to the deployed identification method. The subspace identification methods produce only point estimates but no confidence bounds. The objective for uncertainty quantification is to obtain the parameter estimates and estimates of their covariance from the same dataset. While the statistical properties of estimates from subspace methods have been analyzed in great detail in the automatic control literature in the past, e.g. in [29]–[31] regarding the asymptotic normality and theoretical covariance expressions, these expressions cannot be directly used for an actual covariance estimation in practical applications, since they require in addition e.g. the estimation of the unknown states and their covariances, which are not computed in the estimation of (A, C) or (B, D). A different approach was proposed in [32], where the covariance of estimated parameters is computed easily from the sample covariances of the underlying data covariances and their related sensitivities. The schematics of the framework are presented in Fig. 2.1.

The sensitivity-based covariance propagation is a simple and powerful tool for uncertainty quantification, and is theoretically justified by the statistical delta method [33]. It states that a function of an asymptotically Gaussian variable is also asymptotically Gaussian if its sensitivity is non-zero, and gives the respective covariance expression. Since the data
<figure><img src="image_9.png" /><figcaption>Figure 2.1 – Framework of subspace identification and uncertainty quantification</figcaption></figure>

covariances that are the basis of any subspace method are asymptotically Gaussian, this strategy can be used for characterizing the statistical distributions and in particular the covariance of the state-space matrices, the associated modal parameters and functions thereof, as outlined in the following. Let $\hat{\mathcal{R}}$ be a vector containing all data covariance estimates involved in the chosen subspace method (based on $N$ data samples), then the Central Limit Theorem (CLT)

$$ \sqrt{N}(\hat{\mathcal{R}} - \mathcal{R}) \rightarrow \mathcal{N}(0, \Sigma_{\mathcal{R}}) \qquad (2.1) $$

holds, where an estimate $\hat{\Sigma}_{\mathcal{R}}$ of the covariance can be easily evaluated by the sample covariance based on partitions of the available data. The delta method states that a function $\hat{Y} = f(\hat{\mathcal{R}})$ is also asymptotically Gaussian with

$$ \sqrt{N}(\hat{Y} - Y) \rightarrow \mathcal{N}(0, \mathcal{J}_{Y,\mathcal{R}}\Sigma_{\mathcal{R}}\mathcal{J}_{Y,\mathcal{R}}^T), \qquad (2.2) $$

where $\mathcal{J}_{Y,\mathcal{R}}$ is the derivative of the function with respect to $\mathcal{R}$. This derivative can be obtained from perturbation theory: For a first-order perturbation it holds $\Delta Y \approx \mathcal{J}_{Y,\mathcal{R}} \Delta \mathcal{R}$. Hence, perturbing the functional relationship between $\mathcal{R}$ and $Y$ analytically and neglecting higher-order terms yields the desired derivative, in particular for cases where the functional relationship is not explicit like for the SVD or eigenvalue decomposition. Subsequently, covariance expressions for the estimates satisfy

$$ \hat{\Sigma}_{Y} \approx \hat{\mathcal{J}}_{Y,\mathcal{R}} \hat{\Sigma}_{\mathcal{R}} \hat{\mathcal{J}}_{Y,\mathcal{R}}^T. \qquad (2.3) $$

A bottleneck in such a covariance computation [32] is a numerical one: to evaluate the covariance like $\hat{\Sigma}_Y$ of the involved matrices in the subspace identification methods, the covariance of each matrix entry needs to be evaluated, leading to covariance matrices whose dimensions are squared with respect to the underlying matrices. This leads to a size explosion and an extensive computational burden already for moderate model orders and numbers of sensors, while for realistic problems entire stabilization diagrams at many successive and high model orders need to be evaluated. Furthermore, at the beginning of this line of work, only covariance expressions for the output-only covariance-driven subspace algorithm were available.
## Contributions

My contributions are related to numerical efficiency of uncertainty quantification for subspace methods, which enabled their application on realistic problems for the first time and opened up the development of uncertainty quantification for various subspace methods and related applications in vibration analysis:

* Numerically efficient reformulation of the uncertainty quantification for the covariance-driven subspace method from [32]. The algorithm exploits the structure of the initial sample covariance estimate, optimizes the numerical operations and efficiently exploits the structure of the identification algorithm at multiple model orders, which makes the computations significantly faster and feasible for memory and realistic problem sizes (Section 2.2).

Context: Collaboration with Laurent Mevel.

* Covariance expressions for data-driven and input/output subspace methods for matrices (A, C) for modal parameter estimation. This is the first considerable methodological extension of the uncertainty quantification framework (Section 2.3.1).

Context: Collaboration with Philippe Mellinger during his PhD at Inria [34], and with Laurent Mevel.

* Covariance estimation for matrices (B, D), when some of the input measurements are available, and for the parametric transfer function (Section 2.3.2).

Context: Postdoc of Szymon Greś at Inria, together with Laurent Mevel.

* Uncertainty quantification of modal indicators like the MAC or MPC. We have shown that the preceding first-order frameworks are insufficient for uncertainty quantification for this case, since the indicators are bounded in the interval [0, 1]. To quantify their uncertainties, we have developed an innovative second-order framework (Section 2.4).

Context: Collaboration with Szymon Greś during his PhD at Aalborg University, and with Laurent Mevel.

## 2.2 Efficient multi-order computation

The algorithm for uncertainty quantification of modal parameters from covariance-driven subspace identification has been proposed in [32], where the involved analytical sensitivities for the covariance propagation as well as the initial sample covariance are derived. However, in its direct implementation the size of the involved covariance matrices is considerable, which makes it computationally taxing and causes memory problems even for moderately sized problems. To alleviate this problem, we have made a mathematical reformulation of the algorithm that leads to a memory efficient and fast computation scheme. The key points are:

* The structure of the initial sample covariance estimate $\hat{\Sigma}_R$ (see Eq. (2.3)) is exploited. In practice, the number of data blocks for its estimation is often lower than its dimension, leading to the decomposition $\hat{\Sigma}_R = TT^T$ with an often much smaller matrix $T$.
* The computation of the sensitivities $J_{Y,R}$ is reformulated, avoiding Kronecker products with full matrices (that appear due to vectorization of the involved matrices), but instead exploiting sparse structures. Only products of type $J_{Y,R}T$ are computed, which involve the small matrix $T$, but the full covariance in Eq. (2.3) is never explicitly computed until the final step.

* The covariance computation is optimized for multiple model orders in the stabilization diagram, where the fact is exploited that the first columns of the observability matrix estimate at a higher model order are identical to the observability matrix estimate at a lower model order due to the SVD in (1.4). This has been made explicit for the pseudo-inverse based computation of $A$ in [J5] and for the QR-based computation of $A$ based on Theorem 1.1 in [CA5]. In this way, redundant operations in the computation of the stabilization diagram uncertainties can be avoided.

With this strategy, the algorithmic speed increased by two orders of magnitude of the maximal model order from $O(n_{max}^6)$ to $O(n_{max}^4)$, resulting in computation times of less than a minute for typical problem sizes. Computation times are illustrated in Fig. 2.2, where the conventional implementation fails at a low model order being out of memory, while the efficient implementation scales well and is considerably faster. The applied strategies thus enable the implementation of uncertainty quantification for realistic problem sizes. They will also form the basis for further developments on uncertainty quantification in the next sections as well as in Chapter 4.

In Fig. 2.3, the stabilization diagram from a dataset of the Z24 Bridge is shown, where the estimated standard deviations of the frequencies are shown as horizontal bars. Putting a threshold on the coefficient of variation of the frequencies (standard deviation divided by the frequency) cleans the diagram considerably, as seen in Fig. 2.3 (right).

The uncertainty information of the modal alignments in the diagram can then be used to obtain global mode estimates by a weighted average of the modes in the alignments, where the weights are the inverse covariances [CS48], and to actually obtain the modal alignments

<figure><img src="image_10.png" /><figcaption>Figure 2.2 – Computation times for covariance computation of modal parameters and their point estimates only for different maximal model orders $n_m$. Algorithm 3 denotes the direct implementation, Algorithm 4 denotes the efficient implementation.</figcaption></figure>
<figure><img src="image_11.png" /><figcaption>Figure 2.3 – Stabilization diagram of Z24 Bridge without (left) and with threshold on standard deviation of frequencies (right).</figcaption></figure>

by automated statistical clustering approaches [CS74]. The identified modal parameters and their uncertainties are an important input for damage diagnosis methods, see also Sections 3.2 and 3.3.

[J5] M. Döhler and L. Mevel, “Efficient multi-order uncertainty computation for stochastic subspace identification,” *Mechanical Systems and Signal Processing*, vol. 38, no. 2, 346–366, 2013

[CA5] M. Döhler, X.-B. Lam, and L. Mevel, “Multi-order covariance computation for estimates in stochastic subspace identification using QR decompositions,” in *19th IFAC World Congress, Cape Town, South Africa, 2014*

[CS20] M. Döhler, P. Andersen, and L. Mevel, “Operational modal analysis using a fast stochastic subspace identification method,” in *Proc. 30th International Modal Analysis Conference (IMAC)*, Jacksonville, FL, USA, 2012

[CS48] M. Döhler, P. Andersen, and L. Mevel, “Variance computation of modal parameter estimates from UPC subspace identification,” in *Proc. 7th International Operational Modal Analysis Conference (IOMAC)*, Ingolstadt, Germany, 2017

[CS76] M. Döhler, “Quantification of statistical uncertainties in subspace-based operational modal analysis and their applications,” in *Proc. 9th International Operational Modal Analysis Conference (IOMAC)*, Vancouver, Canada, 2022

## 2.3 Uncertainty in combined deterministic-stochastic identification

The family of subspace methods is big, including covariance-driven and data-driven methods, methods for output-only (stochastic) system identification, and methods for input/output (combined deterministic-stochastic) identification for the case where some of the inputs are
known [1]. They differ in the construction of matrix $\hat{\mathcal{H}}$ from output (and, if available, input) data Hankel matrices such that the observability matrix defines a basis of its column space, see Section 1.2. Different ways of computing $\hat{\mathcal{H}}$ hence lead to different asymptotic covariance estimates, for which the above uncertainty quantification strategy based on the covariance-driven output-only subspace method is generalized. The estimation of $(A, C)$ and $(B, D)$ can be done in two separate steps, and the uncertainty quantification follows this scheme.

### 2.3.1 Matrices A and C

The main difficulty in the extension to other subspace methods lies in the link between the matrix $\hat{\mathcal{H}}$ and the data-related covariance matrices that are used in its construction, in order to perform the first step of uncertainty propagation.

The covariance-driven subspace method is the simplest of the subspace methods for uncertainty quantification, since it only depends on one data-related covariance matrix $Y^+Y^{-T}$, where matrix $Y^-$ is defined analogously to (1.19) but with a past time horizon. This product is a consistent estimate that converges under stationarity assumptions, and therefore allows uncertainty quantification based on CLT (2.1), whereas the dimensions of the data Hankel matrices $Y^-$ and $Y^+$ depend on the number of samples $N$ and do not possess any convergence properties. Herein lies the challenge when going to (output-only or input/output) ‘data-driven’ methods like UPC [1], [7] or N4SID [1], [9], where the matrix $\hat{\mathcal{H}}$ is not immediately linked to output covariances, but is a matrix that grows with the number of data samples. However, it is easy to show that the observability matrix can be equivalently related to $\hat{\mathcal{H}}\hat{\mathcal{H}}^T$ for the purpose of uncertainty propagation, which is again ‘covariance-driven’ and usable for a consistent computation of the covariances. For example, for the case of UPC where $\hat{\mathcal{H}} = Y^+Y^{-T}(Y^{-}Y^{-})^{-1}Y^{-}$, the matrix $\hat{\mathcal{H}}\hat{\mathcal{H}}^T$ solely depends on two data-related covariance matrices, namely $Y^+Y^{-T}$ and $Y^{-}Y^{-}$. In the presence of known inputs, matrix $\hat{\mathcal{H}}$ depends furthermore on covariance matrices between outputs and inputs, and between inputs. To allow a consistent computation of these matrices, the requirement of the inputs to be a quasi-stationary sequence needs to be made, as e.g. in [35]. Under these circumstances, the sample covariance of the involved data covariance matrices can be evaluated and propagated to $\hat{\mathcal{H}}$ (for covariance-driven methods) or $\hat{\mathcal{H}}\hat{\mathcal{H}}^T$ (for data-driven methods). Once this is done, the uncertainty propagation to the observability matrix, matrices $(A, C)$ and the modal parameters is identical as described in Section 2.2, and uses the same strategies for an efficient computation as described therein.

We have illustrated the related uncertainty quantification schemes for several subspace methods in [J12], namely

* Output-only orthogonal projection data-driven algorithm (UPC) [1], [7],

* Input/output covariance-driven algorithm [36],

* Input/output orthogonal projection data-driven (similar to MOESP) [1],

* Input/output oblique projection data-driven (N4SID) [1], [9].

The developed strategies also allow an easy extension to other subspace methods.
### 2.3.2 Matrices *B* and *D*

With the availability of input data, the system matrices (*B*, *D*) of the related input/output state space model (1.3) can be identified in addition to matrices (*A*, *C*) with the method developed in Section 1.5. This allows, amongst others, the identification of the parametric transfer function of the system.

The uncertainty quantification of (*B*, *D*) with a practical and efficient algorithm similar to those developed for (*A*, *C*) was missing in this context. Since estimates of (*A*, *C*) are used in the estimation of (*B*, *D*) in Section 1.5, the uncertainty quantification depends on the output and input/output covariance matrices that are used in the chosen subspace algorithm for the identification of (*A*, *C*), as well as the covariance matrices $\hat{R}_1$ and $\hat{R}_2$ that are used in the estimation of (*B*, *D*) in Section 1.5. All these contributions are carefully gathered for the derivation of the associated sensitivities in the technical development of the algorithm, propagating the uncertainties to (*B*, *D*), and to the phase and magnitude of the transfer function. Computed confidence intervals of the transfer function in comparison to empirical Monte Carlo results are illustrated in Fig. 2.4.

<figure><img src="image_12.png" /><figcaption>Figure 2.4 – First component of phase and magnitude of transfer function $H(z)$ with Monte Carlo and delta method-based confidence intervals.</figcaption></figure>

[J12] P. Mellinger, M. Döhler, and L. Mevel, “Variance estimation of modal parameters from output-only and input/output subspace-based system identification,” *Journal of Sound and Vibration*, vol. 379, 1–27, 2016

[J25] S. Greś, M. Döhler, N.-J. Jacobsen, and L. Mevel, “Uncertainty quantification of input matrices and transfer function in input/output subspace system identification,” *Mechanical Systems and Signal Processing*, vol. 167, 108581, 2022

[CA14] S. Greś, M. Döhler, and L. Mevel, “Variance computation for system matrices and transfer function from input/output subspace system identification,” in *Proc. 21st IFAC World Congress*, 2020
## 2.4 Uncertainty of modal indicators: second-order framework

Modal indicators like the Modal Assurance Criterion (MAC) [37] and the Modal Phase Collinearity (MPC) [38] are widely used in monitoring for the analysis of mode shapes, e.g., for data-model-matching, modal tracking or complexity analysis. They are bounded in the interval [0, 1]. The theoretical value of one is of particular interest, e.g., to infer if the MAC is evaluated for identical mode shapes or if the MPC is evaluated for a real-valued mode shape. Since estimates of the modal indicators are computed from mode shape estimates, they inherit their statistical uncertainties and will never be exactly one, but only close to one. Hence, the quantification of the statistical uncertainties is required in order to evaluate if the modal indicators are sufficiently close to one or not. The statistical characterization of MAC or MPC estimates has never been made before.

Due to their boundedness the asymptotic distribution properties of the modal indicators are very different from the Gaussian distribution of the previously obtained system matrices of the state-space model or the modal parameters. In fact, the first theoretical difficulty for the uncertainty quantification of the modal indicators was the derivation of their distribution properties, which have been unknown before. While the existing framework was based on the (first-order) sensitivity-based propagation of a sample covariance in a Gaussian framework like in Eq. (2.2), we have shown that the first derivative is actually zero for the considered problem class, making this impossible.

The developed second-order framework is novel to the domain and applicable to an entire problem class for which no uncertainty quantification was available before. Its originality is the asymptotic distribution analysis based on a second-order Taylor expansion. For example, when the MAC is evaluated on two mode shape estimates $\hat{\varphi}$ and $\hat{\psi}$, with their theoretical (converged) counterparts being $\varphi_*$ and $\psi_*$ that are equal, then we have shown that the first derivative $\mathcal{J}_{\varphi_*, \psi_*}^{\text{MAC}}$ of the MAC with respect to the mode shapes is zero and the second-order Taylor expansion yields

$$ \text{MAC}(\hat{\varphi}, \hat{\psi}) \approx \underbrace{\text{MAC}(\varphi_*, \psi_*)}_{=1} + \underbrace{\mathcal{J}_{\varphi_*, \psi_*}^{\text{MAC}}(\hat{X} - X)}_{=0} + \frac{1}{2}(\hat{X} - X)^T H_{\varphi_*, \psi_*}^{\text{MAC}} (\hat{X} - X), \quad (2.4) $$

where $H_{\varphi_*, \psi_*}^{\text{MAC}} \in \mathbb{R}^{4r \times 4r}$ is the Hessian, i.e. the second derivative of MAC($\varphi, \psi$) in $\varphi_*$ and $\psi_*$, and $\hat{X}$ is a vector containing the real and imaginary parts of the mode shape estimates that are computed on $N$ data samples and that are asymptotically Gaussian with $\hat{X}_N = \sqrt{N}(\hat{X} - X) \to \mathcal{N}(0, \Sigma_{\varphi_*, \psi_*})$, cf. Eq. (2.1)-(2.3). The Hessian is developed in detail in [J20]. The asymptotic properties of MAC($\hat{\varphi}, \hat{\psi}$) follow from (2.4) as

$$ N(1 - \text{MAC}(\hat{\varphi}, \hat{\psi})) \approx \hat{X}_N^T \bar{H}_{\varphi_*, \psi_*}^{\text{MAC}} \hat{X}_N = Q(\hat{X}_N), \quad (2.5) $$

where $\bar{H}_{\varphi_*, \psi_*}^{\text{MAC}} = -\frac{1}{2}H_{\varphi_*, \psi_*}^{\text{MAC}}$. This equation describes a quadratic form of a Gaussian variable that is directly linked to the MAC estimate. The exact distribution of such a quadratic form is difficult to characterize; it is clearly non-Gaussian, and thus cannot be defined solely by its mean and variance. However, an approximation of the distribution using a scaled and shifted $\chi^2$ distribution [39] can be applied when the matrix $\bar{H}_{\varphi_*, \psi_*}^{\text{MAC}}$ is positive semi-definite,
which we have shown to be true. With this approximation, the statistical properties of the quadratic form can be related to the asymptotic covariance $\Sigma_{\varphi_*, \psi_*}$ of the mode shape estimates, and subsequently the statistical properties of the MAC estimate can be quantified through relation (2.5).

This leads to a practical approximation of the MAC distribution by a scaled and shifted $\chi^2$ distribution whose confidence bounds can be easily computed. To illustrate the distribution properties, Fig. 2.5 shows histograms of the MAC and the developed distribution fits when evaluated on mode shape estimates of different (left) as well as on equal modes (right). The resulting statistical evaluation of modal indicators can be highly relevant for many applications in vibration analysis, e.g., to evaluate if a change is related to damage, or in model updating.

<figure><img src="image_13.png" /><figcaption>Figure 2.5 – Distribution fits of MAC estimates from delta method, together with histograms from Monte Carlo simulation. Left: MAC of estimates of different mode shapes with Gaussian approximation, right: MAC of estimates of equal mode shape with scaled and shifted $\chi^2$ approximation.</figcaption></figure>

[J20] S. Greś, M. Döhler, P. Andersen, and L. Mevel, “Uncertainty quantification for the modal phase collinearity of complex mode shapes,” *Mechanical Systems and Signal Processing*, vol. 152, 2017, pp. 107436–107451.

[J21] S. Greś, M. Döhler, and L. Mevel, “Uncertainty quantification of the modal assurance criterion in operational modal analysis,” *Mechanical Systems and Signal Processing*, vol. 152, 2017, pp. 107457–107471.

[CS58] S. Greś, M. Döhler, P. Andersen, and L. Mevel, “Variance computation of the modal assurance criterion,” in *Proc. 28th Conference on Noise and Vibration Engineering (ISMA)*, Leuven, Belgium, 2018.

[CS60] S. Greś, M. Döhler, P. Andersen, and L. Mevel, “Variance computation of MAC and MPC for real-valued mode shapes from the stabilization diagram,” in *Proc. 8th International Operational Modal Analysis Conference (IOMAC)*, Copenhagen, Denmark, 2019.
## 2.5 Conclusions and outlook

With these contributions, uncertainty quantification for subspace-based system identification has become feasible for realistic problem sizes in practice, and thus enabled further methodological development for the whole family of subspace methods, including data- and covariance-driven techniques as well as output-only (stochastic) and input/output (combined deterministic-stochastic) techniques. The subsequent uncertainty quantification of the derived modal indicators is particularly original as it employs a second-order framework (and non-Gaussian distributions) for the first time in asymptotic vibration analysis.

The applicability of the developed methods to real problems is also shown by technology transfers to commercial software. The efficient uncertainty computation (Section 2.2) adapted for the data-driven UPC subspace method (Section 2.3.1) and the methods for modal indicator uncertainty quantification (Section 2.4) have been transferred to ARTEMIS Modal [28]. Furthermore, the efficient uncertainty quantification for covariance-driven subspace identification has been transferred to the French seismic data acquisition system provider SERCEL with the purpose of developing an integrated SHM solution with their high-end sensors.

These activities are ongoing, and we explore the uncertainty quantification within robust and automated modal parameter estimation under environmental changes in an industrial ANR project [CS74] with postdoc Johann Priou. With Szymon Greś, we went back to the roots of the uncertainty quantification methods to consider the effects of unknown and misspecified model order in the uncertainty computations, which turns out to be an important issue as the sensitivities of image and null space of a matrix are different [CS79], [CS82]. Further ongoing work concerns the uncertainty quantification of FRF's and input/output system identification in the frequency domain, in collaboration with Hottinger Brüel & Kjær's PhD student Mikkel Steffensen.
# Uncertainty quantification for deterministic damage diagnosis

## 3.1 Introduction

Quantification of statistical uncertainty in damage diagnosis is essential for statistical decision-making. Amongst others, it also helps to assess the equivalence between parameters of a numerical model and the related estimated parameters obtained from data. For damage diagnosis, data from a current state is typically compared to a data-driven and/or model-based reference. Typical steps in the process include the operational evaluation of the structure, the acquisition of system response data through permanently installed sensors, the extraction of damage-sensitive features from the data, and their subsequent statistical evaluation [40]. The diagnosis levels for the evaluation are classically divided into detection, localization, quantification, and lifetime prognosis.

Methods for damage detection are the most developed since they can operate purely data-based and do not require a physical model of the monitored structure. Examples are algorithms from pattern classification and statistical process control, which usually consider the statistical uncertainty related to the data [40]. However, in general no physical interpretation is possible when operating purely data-driven. Methods for damage localization and the higher levels are more sophisticated since some link between the measurement data and the physical properties of the structure is required in order to localize and quantify the physical change in the structure due to damage. Such a link is usually given by a finite element (FE) model, which can be used to link the damage-sensitive feature from the data with a model-based parameterization for a physical interpretation.

Typically, methods for damage localization or quantification use the extracted damage-sensitive feature (like the modal parameters) in a deterministic way, and rarely consider its statistical uncertainty when inferring about the damage location or extent [41]–[43]. The
methods are often solely based on mechanical or physical grounds. For example, modal pa-
rameters are computed and used for FE model updating to infer about the changed param-
eters and their extent, but without consideration of their covariance. Often, these problems
are ill-posed or ill-conditioned, since different parameter changes may yield a similar change
in the feature, which is even aggravated by the presence of uncertainties in the feature. These
uncertainties are unavoidable due to unknown excitation and measurement noise, and can be
evaluated from system identification, as shown in Chapter 2. Hence, it seems natural to also
use the available information on the uncertainties of the used estimates, with the goal of en-
hancing the damage diagnosis methods. In contrast to Bayesian methods that assume a prior
distribution of the model parameters and sometimes also of the modal parameter estimates
in this context [43], distribution assumptions are not necessary in our work, and the actual
modal parameter uncertainty is propagated as it is obtained from the data. In this chapter
two examples are given where the quantified modal parameter uncertainties are integrated
into originally deterministic methods, namely for damage localization and for finite element
model updating.

The considered damage localization method is the *damage locating vector* (DLV) approach
[44], [45], where a vector is estimated in the null space of the transfer matrix difference be-
tween reference and damaged states, which is computed from the modal parameter estimates.
Theory shows that when applying this vector as a virtual load to the FE model of the struc-
ture, the resulting stress over the damaged elements is zero. Such an approach is particularly
promising, as the problems of damage localization and quantification are divided into sepa-
rate steps, which alleviates the typical ill-posedness to some extent. However, the computed
stresses that serve as a damage location indicator are deterministic, and we have shown that
the integration of the estimation uncertainties into the method and the subsequent statisti-
cal evaluation increases the performance of the approach significantly. Regarding FE model
updating, we have considered objective functions to match modal parameter estimates to
their model counterparts. Then, the uncertainty of the updated model parameters due to
the modal parameter uncertainty can be evaluated by uncertainty propagation. Furthermore,
the uncertainties of the modal parameters can be directly considered in the formulation of
the objective function as well as in the stopping criterion for the optimization, resulting in
an efficient method that can be used for damage quantification.

Contributions

My contributions lie in the uncertainty quantification of quantities that are related to damage
localization and quantification, based on the uncertainty of modal parameter estimates, and
the subsequent statistical decision making for damage diagnosis:

* Integration of the estimation uncertainties of the modal parameters into the DLV framework, in order to develop a robust statistical evaluation for damage localization. In a second step, further constraints of the basic method were alleviated through statistical decision making to improve its applicability, including a follow-up considering temperature robustness (Section 3.2).

Context: My postdoc at Northeastern University, Boston, with Dionisio Bernal; PhD of Luciano Marin [46] (supervision together with Laurent Mevel and Dionisio Bernal);
PhD of Delwar Bhuyan [47] (supervision together with Franck Schoefs and Yann Lecieux at University of Nantes, and Laurent Mevel); Delwar Bhuyan's postdoc at BAM with Falk Hille and in collaboration with Jean Dumoulin at Université Gustave Eiffel together with postdocs Nicolas Le Touz and Guillaume Gautier.

* Statistical evaluation of updated model parameters in FE updating, based on the modal parameter uncertainties, and integration of modal parameter uncertainties in the optimization procedure to conceive a more robust and efficient updating method for damage quantification (Section 3.3).

Context: Postdoc of Guillaume Gautier at Inria; collaboration with Roger Serra and Jean-Mathieu Mencik at INSA Val de Loire, and Laurent Mevel; postdoc of Szymon Greś at Inria.

## 3.2 Statistical damage localization with Damage Locating Vector approach

The damage locating vector (DLV) approach systematically interrogates changes in the structural flexibility in order to locate damages [44]. In its generalization to output-only data [45], [48], an analogous damage locating vector is estimated in the null space of the transfer matrix difference $\delta G(s)$ between reference and damaged states, which can be entirely obtained from the modal parameter estimates in the reference and in the damaged states. However, the transfer matrix

$$G(s) = C(sI - A_c)^{-1}B_c + D_c \in \mathbb{C}^{r \times r},$$

where $A_c$, $B_c$ and $D_c$ denote the system matrices of the continuous-time system, cannot be estimated from output-only measurements. It holds [45], [48]

$$G(s) = R(s)D_c, \quad \text{where} \quad R(s) = C_c(sI - A_c)^{-1} \begin{bmatrix} C_c A_c \\ C_c \end{bmatrix}^{\dagger} \begin{bmatrix} I \\ 0 \end{bmatrix}. \qquad (3.1)$$

under the condition that the system order satisfies $n \le 2r$, i.e., the number of identified modes is bounded by the number of sensors. The difference between the transfer matrices in both damaged (variables with tilde) and healthy states is $\delta G(s) = \tilde{G}(s) - G(s)$. Assume that damage is due to changes in stiffness and mass is constant. Then $D_c = \tilde{D}_c$, and the matrices $\delta G(s)$ and $\delta R(s)^T = \tilde{R}(s)^T - R(s)^T$ have the same null space [45]. The desired load vector $v(s)$ is obtained from the null space of the $\delta R(s)^T$ with an SVD. To compute the stress field over the elements of the structure, the load vector $v(s)$ is applied to the FE model. First, the vector $v(s)$ is expanded to load vector $\nu(s)$ at all DOFs of the model, whose entries are those of $v(s)$ at the sensor coordinates and zero elsewhere. From this vector, the nodal displacements $y(s) = (Ms^2 + Cs + K)^{-1}\nu(s)$ are computed based on the FE model at all DOFs, from which stress resultants are evaluated for each structural element of the model and stacked into stress vector $S(s)$. This relation between stress and load is linear and can be expressed by a matrix $L_{\text{model}}(s) \in \mathbb{C}^{l \times r}$ based on the FE model of the structure, where $l$ is the number of considered elements, satisfying

$$S(s) = L_{\text{model}}(s)v(s). \qquad (3.2)$$
The stress vector $S(s) \in \mathbb{C}^l$ indicates potential damage for elements with corresponding entries in $S(s)$ that are close to zero [44], [45]. However, due to modal truncation errors and the estimation uncertainty, the stress field over damaged elements is never exactly zero. However, in the original formulation of the method, the stress vector was evaluated deterministically. Instead, we have proposed a statistical evaluation, which is based on the modal parameter uncertainties as obtained in the previous chapter.

In a first step, the modal parameter uncertainties are propagated to matrices $A_c$ and $C$ in the modal basis, and then to $R(s)$ in (3.1) in both states. After developing the null space uncertainties for the SVD, the propagation continues to $v(s)$ and finally to $S(s)$ in (3.2). Therein, the extraction of matrix $\mathcal{L}_{\text{model}}(s)$ from an FE model may not be straightforward, but can always be achieved by computing stresses from unit loads at the sensor DOFs, and assembling them as columns of $\mathcal{L}_{\text{model}}(s)$. The efficient numerical strategies from the last chapter are adopted in the uncertainty computation, and lead to covariance $\Sigma_t = \text{cov}(S_t)$ of the stress(es) corresponding to an element $t$. The statistical evaluation of the stress is finally made in the test statistic

$$ \chi_t^2 = S_t^T \Sigma_t^{-1} S_t \qquad (3.3) $$

for each structural element $t$. An example is given in Fig. 3.1 for damage localization on a beam in the lab, where holes were drilled in the oval zone indicated in the figures. The deterministic stresses and the statistical counterparts are visualized in the color maps in the figures, where red corresponds to the smallest values. It can be seen that successful localization became possible thanks to the statistical evaluation.

The formulation (3.3) offers a natural weighting of different stress components by their uncertainties for a meaningful evaluation. Indeed, several stress resultants can be computed for an element depending on its nature (plate element, 2D beam, 3D beam, ...). Moreover, the stresses can (and should) be computed for different values of the parameter $s$, which is recommended to be in the vicinity of the identified modes [45]. The stresses can also be

<figure><img src="image_14.png" /><figcaption>Figure 3.1 – Localization results of beam experiment: joint stress computation and statistical evaluation at ten s-values and different damage levels. The oval indicates the actual damage location.</figcaption></figure>
computed based on different mode sets that satisfy the condition that there are not more
modes than sensors. Thus, the statistical evaluation is the key for application of the method,
in particular when the choice of those parameters is more or less automated, as the provided
weighting will naturally emphasize values with low uncertainty and exclude values with high
uncertainty when merging the different computed stresses for each structural component.

The method has also been extended the case of changing temperature, where the tempera-
ture influence is removed from the modal parameter estimates based on sensitivities obtained
from a temperature-dependent parameterization of the FE model. Here it also has been
shown that the performance of the damage localization is largely increased when considering
the modal parameter uncertainties, as illustrated in Fig. 3.2. It can be seen that the rate of
successful damage localization is always higher when the considering the uncertainties in the
statistical tests compared to the deterministic evaluation of the stress estimates.

<figure><img src="image_15.png" /><figcaption>Figure 3.2 – Success rates for DLV-based damage localization with and without considering modal parameter uncertainties for different temperature scenarios. Success rates are higher when considering the uncertainties in the evaluation (red vs. black curve, pink vs. blue curve).</figcaption></figure>

[J6] M. Döhler, L. Marin, D. Bernal, and L. Mevel, “Statistical decision making for damage localization with stochastic load vectors,” *Mechanical Systems and Signal Processing*, vol. 39, no. 1-2, 426–440, 2013

[J11] L. Marin, M. Döhler, D. Bernal, and L. Mevel, “Robust statistical damage local-
ization with stochastic load vectors,” *Structural Control and Health Monitoring*,
vol. 22, no. 3, 557–573, 2015

[J15] M. Bhuyan, M. Döhler, Y. Lecieux, L. Mevel, and F. Schoefs, “Statistical damage localization with stochastic load vectors using multiple mode sets,” *Structural Health Monitoring*, vol. 16, no. 5, 518–535, 2017
[J17] M. Bhuyan, G. Gautier, N. Le Touz, M. Döhler, F. Hille, J. Dumoulin, and L. Mevel, "Vibration-based damage localization with load vectors under temperature changes," Structural Control and Health Monitoring, vol. 26, no. 11, e2439, 2019

[CA4] L. Marin, M. Döhler, D. Bernal, and L. Mevel, "Uncertainty quantification for stochastic damage localization for mechanical systems," in Proc. 8th IFAC Symposium on Fault Detection, Diagnosis and Safety of Technical Processes (SAFE-PROCESS), Mexico City, Mexico, 2012

[CA10] M. Bhuyan, G. Gautier, M. Döhler, Y. Lecieux, L. Mevel, and F. Schoefs, "Statistical damage localization in mechanical systems based on load vectors," in Proc. 20th IFAC World Congress, Toulouse, France, 2017

[CS26] L. Marin, M. Döhler, D. Bernal, and L. Mevel, "Damage localization using a statistical test on residuals from the SDDLV approach," in Proc. 31st International Modal Analysis Conference (IMAC), Garden Grove, CA, USA, 2013

[CS44] M. Bhuyan, M. Döhler, and L. Mevel, "Statistical damage localization with stochastic load vectors using multiple mode sets," in Proc. 8th European Workshop on Structural Health Monitoring (EWSHM), Bilbao, Spain, 2016

[CS61] M. D. H. Bhuyan, N. Le Touz, G. Gautier, M. Döhler, F. Hille, J. Dumoulin, and L. Mevel, "Load vector based damage localization with rejection of the temperature effect," in Proc. 8th International Operational Modal Analysis Conference (IOMAC), Copenhagen, Denmark, 2019

## 3.3 Uncertainties in finite element model updating

Modal parameter estimates are often used for FE model updating in order to adjust the model parameters such that modal parameters of the model match the actual estimates [43]. One aim of FE model updating is to locate and quantify damage based on the changed model parameters.

Since the updated model parameters rely on the modal parameter estimates, it can be important to evaluate the uncertainty of the obtained model parameters that is due to the estimation uncertainty. We have made such an analysis for model updating with a subspace fitting approach, where the objective function is based on matching the observability matrix estimate from the data with the one from the model,

$$ \hat{\theta} = \operatorname*{argmin} \|\operatorname{vec}(O(\theta) - \hat{O}T)\|^2, \qquad (3.4) $$

where $T = \hat{O}^\dagger O(\theta)$ is the appropriate change of basis matrix. This can be solved iteratively with the Gauss-Newton method. For uncertainty quantification, the uncertainties of the observability matrix estimate (based on the modal parameter uncertainties) can be propagated to the updated parameter vectors in these iterations, with the technical details given in [J14].
The iterative uncertainty propagation was the original part in this work, and we could empirically show their convergence. We have validated the so obtained model parameter covariance in numerical Monte Carlo simulations and on experimental data of a beam.

In the next step, the uncertainties are not only analyzed, but also used in the updating procedure itself. The updating problem is in general ill-posed (or at least ill-conditioned) and small changes in the modal parameter estimates may yield large changes in the model parameters, which is even aggravated by the unavoidable presence of the modal parameter uncertainties. We have proposed to consider the uncertainties directly in the formulation of the objective function as well as in the stopping criterion for the optimization. A more direct objective function was used for this task, namely a function that explicitly indicates the relative distance between the estimated modal parameters and their model counterparts. It is expressed as the sum of the respective differences $\Delta_{f_i}^{\theta} = 1 - f_i^{\theta}/\hat{f}_i$ and $\Delta_{\text{MAC}_i}^{\theta} = 1 - \text{MAC}(\varphi_i^{\theta}, \hat{\varphi}_i)$ for all considered modes $i = 1, \dots, m$, as

$$F(\theta) = \sum_{i=1}^{m} |\Delta_{f_i}^{\theta}| + \sum_{i=1}^{m} \Delta_{\text{MAC}_i}^{\theta}. \qquad (3.5)$$

Including the uncertainties, model frequencies and mode shapes (via the MAC) are penalized in the objective function when they are outside the confidence bounds of their estimated counterparts by setting $\Delta_{f_i}^{\theta}$ or $\Delta_{\text{MAC}_i}^{\theta}$ to 1, respectively. For this, the evaluation of the MAC uncertainties is based on the derivations in Section 2.4, which are adapted here to the case where the MAC is computed between a mode shape estimate and a (deterministic) mode shape from a model. This leads to a steeper objective function, as illustrated in Fig. 3.3. To minimize (3.5), the Covariance Matrix Adaptation Evolution Strategy (CMA-ES) [49] is applied, and adopted to consider the confidence bounds of the modal parameter estimates: the optimization search can be terminated once the model-based modal parameters are within the confidence bounds of their data-based counterparts, leading to an efficient algorithm

<figure><img src="image_16.png" /><figcaption>Figure 3.3 – Modal parameter-based objective function to be minimized in the model updating, for two model parameters in the vicinity of their true values. Without (left) and with consideration of modal parameter uncertainties (right).</figcaption></figure>
for the updating problem. When this methodology is used for damage quantification, a localization should be performed beforehand (like from the previous section) to narrow down the possibly changed parameters. Then, quantification can only be carried out with respect to distinguishable parameters, i.e., a clustering procedure should be carried out before that regroups parameters with (nearly) the same influence on the objective function.

In Fig. 3.3, the modal parameter-based objective function is illustrated without and with the consideration of the uncertainties. First of all, it can be seen that the objective function is indeed steeper when considering the uncertainties, which should facilitate the optimization search. Second, it can be seen that the minimum of the objective function (red line) that is computed on the modal parameter estimates is not taken on at the true parameter values (green line) due to the modal parameter uncertainties, but just close to it. This also confirms that convergence of the optimization search to the minimum of the objective function is actually not required, and it is sufficient to stop the search once arrived in the yellow region indicating the values of the objective function where the model-based modal parameters are within the confidence bounds of their estimated counterparts.

[J14] G. Gautier, L. Mevel, J.-M. Mencik, R. Serra, and M. Döhler, "Variance analysis for model updating with a finite element based subspace fitting approach," *Mechanical Systems and Signal Processing*, vol. 91, 142–156, 2017

[J22] S. Greś, M. Döhler, and L. Mevel, "Statistical model-based optimization for damage extent quantification," *Mechanical Systems and Signal Processing*, vol. 160, 107894, 2021

[CS42] G. Gautier, L. Mevel, J.-M. Mencik, M. Döhler, and R. Serra, "Parameter uncertainties quantification for finite element based subspace fitting approaches," in *Proc. 8th European Workshop on Structural Health Monitoring (EWSHM)*, Bilbao, Spain, 2016

[CS70] S. Greś, M. Döhler, and L. Mevel, "Statistical optimization for subspace-based damage quantification," in *10th International Conference on Structural Health Monitoring of Intelligent Infrastructure (SHMII)*, 2021

## 3.4 Conclusions and outlook

This chapter illustrates how the usage of estimation uncertainties in originally deterministic approaches for damage localization or quantification leads to improved methods with better performance. The methods gain robustness towards noise, which increases their applicability on real-world problems.

The application on a sophisticated laboratory test case is ongoing, namely on a mock-up of the Saint-Nazaire Bridge at GeM, Nantes University, in collaboration with Yann Lecieux [CS83]. The goal is to use this experiment as a benchmark for damage diagnosis, where the damages are cable failures of the cable-stayed bridge.
# Sensitivity-based statistical damage diagnosis

## 4.1 Introduction

In this chapter, the damage diagnosis problem is considered in a statistical framework from the outset, where a damage-sensitive feature is computed from measurement data and analyzed by statistical hypothesis testing. The link to physical system properties is made by the feature's sensitivities with respect to a physical parameterization. In contrast, the consideration of statistical uncertainties and physical system properties was made the other way round in Chapter 3, where data-related statistical uncertainties were included in damage diagnosis methods that were already based on physical properties.

Two key points of the methods in this chapter are indeed the intrinsic consideration of statistical uncertainties in the damage-sensitive feature, and the consideration of the model parameters through model-based sensitivities. This allows noise robustness in the feature evaluation (in contrast to purely deterministic methods) on the one hand, and a simple link between data and physical parameters for damage diagnosis (without the need of model updating) on the other hand. The considered sensitivities are valid for small parameter changes (or features that linearize well with respect to the parameters), which is just a small restriction in practice since the damage diagnosis is particularly relevant when damage is incipient.

For detection, the question is if a change in the damage-sensitive feature is statistically significant or not. For localization, the question is which physical parameters are responsible for the change in the feature, which is answered by hypothesis tests with respect to each of these physical parameters, usually related to a finite element (FE) model. This way, the typical ill-conditioning of the inverse problem is avoided, since the number of physical parameters is typically large with respect to the number of sensors/the information content
that can be extracted from measurements. Once the changed parameters are narrowed down, damage quantification can be performed in those parameters only.

The asymptotic local approach framework to change detection [50] offers a statistical framework for the associated feature evaluation and diagnosis with remarkable theoretical properties, the most important one being the characterization of the probability distributions of the feature in both reference and damaged states by means of the central limit theorem (CLT). While statistical damage detection and localization with subspace-based damage features have been developed with this framework in the past [51]–[54], robustness for typical application cases for vibration monitoring was somehow lacking, e.g., when the ambient excitation properties change between measurements, when the reference model itself has uncertainties (as it is usually obtained from data), or when temperature or other environmental variations perturb the damage diagnosis. Furthermore, robustness for the application of the damage localization with experimental data – beyond numerical simulations – was missing in this framework.

With the asymptotic local approach, the probability distribution of the damage-sensitive features is characterized by asymptotic theory with a CLT, which allows a “batch-wise” analysis of datasets. In contrast, the Kalman filter offers inference on the level of the individual data samples, which also allows an analysis within the datasets. However, for *fault detection and isolation* (corresponding to damage detection and localization) in this context, the considered faults are typically *additive* where an additional fault term appears in the state equation [55], but not *multiplicative* as in our case of vibration monitoring where both the state transition matrix *A* and the states $x_k$ change with a damage. The challenge is then to link parametric system changes to the filter equations. For this, we have developed a sensitivity-based Kalman filter-based framework for statistical damage diagnosis in a more exploratory work, assuming again small parameter changes. Here, the innovations computed on a current dataset with the Kalman filter of the reference system are the “damage-sensitive features”. The originality of this work is the parameterization of the innovation bias in terms of the system parameter change, based on a perturbation analysis using model-based sensitivities.

## Contributions

My main contribution to damage diagnosis in the local approach framework lies in developing theory for typical application cases, ensuring robustness for real applications (Section 4.2):

* Robustness of data-driven features towards naturally changing excitation properties (Section 4.2.2).

Context: My PhD, supervised by Laurent Mevel; my postdoc at BAM Federal Institute of Materials Research and Testing, Germany, with Falk Hille; PhD of Szymon Greś at Aalborg University.

* Robustness towards uncertainties in the data-driven reference model (Section 4.2.3), and consideration of changes due to known nuisance parameters like temperature (Section 4.2.4).
Context: PhD of Eva Viefhues [56] at BAM Federal Institute of Materials Research and Testing, Germany (supervision together with Laurent Mevel and Falk Hille).

* Robustness for damage localization considering sensitivity-based clustering approaches, and its first validation on experimental data (Section 4.2.5).

Context: PhD of Saeid Allahdadian at University of British Columbia, Canada (supervision together with Carlos Ventura), in collaboration with Laurent Mevel.

* Numerically stable and efficient computation of the test statistics in the diagnosis framework. This is an important aspect, since the involved covariance matrices can be ill-conditioned when computed from real data, or model-based sensitivity matrices can be ill-conditioned (Section 4.2.6).

Context: Collaboration with Laurent Mevel and Qinghua Zhang.

Finally, we have developed a sensitivity-based Kalman filter-based framework for statistical damage diagnosis (Section 4.3) together with Qinghua Zhang and Laurent Mevel as an alternative to the previous asymptotic local approach framework.

## 4.2 Robustness for local approach methods

### 4.2.1 Introduction

The asymptotic local approach to change detection has been developed at Inria in the 1980's [50] and offers a universal statistical framework for the evaluation of a residual (i.e., a damage-sensitive feature) with respect to a system parameterization. The approach does not require (usually Gaussian) assumptions on the probability distribution of measurement or modeling errors, but instead has the remarkable ability of transforming quite general residuals with unknown probability distributions into a standard Gaussian framework, thanks to a CLT. The local approach is thus a powerful tool for statistical evaluation of residuals. The essential assumption behind the local approach is that the considered faults (i.e., damages) are characterized by small parametric changes.

Let a system be characterized by a parameter vector $\theta$, which relates e.g. to material properties of different components or parts of the structure. The value of $\theta$ is unknown. Assume that measurements $\mathcal{Y}_N = \{y_1, \dots, y_N\}$ of length $N$ collected from the system are the realization of an asymptotically stationary stochastic process depending on the parameter vector $\theta$. Damage corresponds to deviations of the considered system from the nominal behavior characterized by some reference parameter value $\theta_0$ (whose value is known). Small deviations are considered for the diagnosis of incipient faults, assuming the close hypotheses [50]

$$ 
\begin{aligned}
\mathbf{H}_0: \theta = \theta_0 & \quad (\text{reference system}), \\
\mathbf{H}_1: \theta = \theta_0 + \delta/\sqrt{N} & \quad (\text{faulty system}),
\end{aligned} 
\quad (4.1) 
$$

where vector $\delta$ is unknown but fixed. This corresponds to a very general description of a change in system parameter $\theta$. With this framework, very small changes in $\theta$ can be detected if the number of measurements $N$ is large enough, by linking this change to $1/\sqrt{N}$. Confronting the measurements to the reference system, a residual vector $\zeta(\theta_0, \mathcal{Y}_N)$ with mean
<figure><img src="image_17.png" /><figcaption>Figure 4.1 – Theoretical distribution of the test statistic before and after a damage event</figcaption></figure>

zero in the reference state is designed, satisfying some basic conditions [50] like identifiability, differentiability with respect to $\theta$ in a neighborhood of $\theta_0$, and the underlying process satisfying some mild mixing conditions. Then, the main result of the local approach based on hypotheses (4.1) is the CLT

$$ \zeta(\theta_0, \mathcal{Y}_N) \xrightarrow{d} \begin{cases} \mathcal{N}(0, \Sigma) & \text{under } \mathbf{H}_0 \\ \mathcal{N}(\mathcal{J}\delta, \Sigma) & \text{under } \mathbf{H}_1 \end{cases} \quad (4.2) $$

for $N \to \infty$, where $\mathcal{J}$ and $\Sigma$ are the asymptotic sensitivity at $\theta_0$ and covariance of the residual, respectively. Estimates of $\mathcal{J}$ and $\Sigma$ are obtained based on data in the reference state. Note that the modeling of the change in $\theta$ as $\delta/\sqrt{N}$ in hypothesis (4.1) has been introduced to obtain CLT (4.2) with constant mean $\mathcal{J}\delta$ under $\mathbf{H}_1$, resulting from a first-order Taylor expansion of the residual.

Damage detection, i.e., change detection in $\theta$, corresponds then to deciding if $\delta \neq 0$, and damage localization corresponds to deciding which parameters in $\theta$ are changed, that is which components of $\delta$ are non-zero. Based on (4.1)–(4.2), hypotheses tests in a Gaussian setting can be performed. The generalized likelihood ratio (GLR) test for the detection problem leads to the test statistic [51], [52]

$$ t_{\text{global}} = \zeta^T \Sigma^{-1} \mathcal{J} (\mathcal{J}^T \Sigma^{-1} \mathcal{J})^{-1} \mathcal{J}^T \Sigma^{-1} \zeta, \quad (4.3) $$

which follows a $\chi^2$ distribution with $\nu = \text{rank}(\mathcal{J})$ degrees of freedom and the non-centrality parameter

$$ \lambda = \delta^T F \delta, \quad (4.4) $$

see Fig. 4.1, where $F = \mathcal{J}^T \Sigma^{-1} \mathcal{J}$ is the Fisher information matrix. To decide between $\mathbf{H}_0$ and $\mathbf{H}_1$, the test variable $t_{\text{global}}$ is compared to a threshold.

Damage localization corresponds to deciding if a component $\delta_i \neq 0$. So far, it has been assumed that the other components are unchanged for such a test, i.e., $\delta_i = 0$, leading to the test statistic [52], [54]

$$ t_{\text{sens}}^i = \zeta^T \Sigma^{-1} \mathcal{J}_i (\mathcal{J}_i^T \Sigma^{-1} \mathcal{J}_i)^{-1} \mathcal{J}_i^T \Sigma^{-1} \zeta, \quad (4.5) $$

which is called sensitivity test, where $\mathcal{J}_i$ is the $i$-th column of $\mathcal{J}$. The test is computed for each parameter $\theta_i$, to decide individually if there is a change in a particular parameter component
or not. It should be noted that parameters with close sensitivities are not distinguishable in the tests and should be clustered beforehand [54].

A classical residual that has been used in previous works [51]–[54] is the subspace-based residual

$$ \zeta(\theta_0, \mathcal{Y}_N) \stackrel{\text{def}}{=} \sqrt{N} \text{vec}(S^T \hat{\mathcal{H}}), \quad (4.6) $$

where $\hat{\mathcal{H}}$ is the block Hankel matrix containing the output covariances $R_j = E(y_k y_{k-j}^T)$, and $S = S(\theta_0)$ is the left null space of $\mathcal{H}$ from the reference system that can be estimated from $\hat{S} = \hat{U}_2$ in (1.4). Such a residual has the advantage that it can be directly computed and evaluated on measurement data, without complex computations or modal parameter estimation (and the related issues like mode tracking therein). However, due to the direct computation from data it is also more affected by changing operational and environmental conditions.

To ensure applicability of this universal framework under real operation conditions, several robustness issues are considered in the following sections. They are related to the residual definition, its statistical characterization, or the type of considered hypothesis tests used for the evaluation, amongst others.

### 4.2.2 Robustness to change in ambient excitation

Damage-sensitive features $\zeta$ that can be simply and directly computed from measurement data $\mathcal{Y}_N$ are very convenient for an automated damage diagnosis. However, the more directly they are computed from the data, the more they may be influenced by noise or other nuisance. This either needs to be considered in their evaluation, or the residual needs to be designed to be robust towards such influence.

For example, the mean and the covariance of the output covariances in $\hat{\mathcal{H}}$ in (4.6) are directly influenced by the properties of the ambient excitation. Residual (4.6) is easy to compute; however since the ambient excitation properties may change between measurements (e.g., traffic on a bridge is different during the day) this can lead to false alarms or false negatives if not properly accounted for. Residual (4.6) has been designed to check the orthogonality between the left null space $S$ from the reference state and matrix $\hat{\mathcal{H}}$ from the current state, since a change in the dynamic system properties destroys the orthogonality. An analogous orthogonality check can be made between $S$ and the principal singular vectors $\hat{U}_1$ of $\hat{\mathcal{H}}$ as computed in (1.4), with the advantage that the singular vectors have norm one and thus do not depend on ambient excitation levels. The corresponding robust residual writes

$$ \xi(\theta_0, \mathcal{Y}_N) \stackrel{\text{def}}{=} \sqrt{N} \text{vec}(S^T \hat{U}_1), \quad (4.7) $$

and its covariance and sensitivity have been derived in detail in [J7], [J8] for the computation of the statistical damage detection test (4.3). Besides showing the performance on simulated data, we have applied the approach on experimental data of a steel jacket structure in the lab of BAM, where different levels of white noise excitation have been applied by a shaker, and damage has been introduced by loosening bolts of a flange connection. In Fig. 4.2, the values of the test statistic are shown for different excitation levels and an increasing number of loosened bolts. The threshold is empirically determined on the test values in the reference
<figure><img src="image_18.png" /><figcaption>Figure 4.2 – Comparison of damage detection tests on a steel frame structure under different excitation levels. Left: Test with conventional residual (4.6). Right: Test with robust residual (4.7).</figcaption></figure>

state. It can be seen that the test values are much less influenced by the excitation levels when using the robust residual (Fig. 4.2, right) compared to the conventional residual (Fig. 4.2, left), and that smaller damage (starting from three loosened bolts instead of seven) can be detected.

A different residual can be directly constructed on the Hankel matrices in the reference and current states, without the left null space computation, and evaluated with the Mahalanobis distance by

$$d(\hat{\mathcal{H}}_{\text{ref}}, \hat{\mathcal{H}}_{\text{test}}) = \sqrt{\text{vec}(\hat{\mathcal{H}}_{\text{ref}} - \hat{\mathcal{H}}_{\text{test}})^T \Sigma_{\hat{\mathcal{H}}_{\text{ref}}}^{-1} \text{vec}(\hat{\mathcal{H}}_{\text{ref}} - \hat{\mathcal{H}}_{\text{test}})}. \quad (4.8)$$

The purpose of such a distance definition is to provide a simpler and more direct test computation than with the classical subspace residual above. Of course, there are strong parallels between both, such as the dependence on the ambient excitation properties that perturb the distance computation (4.8) when $\hat{\mathcal{H}}_{\text{ref}}$ and $\hat{\mathcal{H}}_{\text{test}}$ are obtained under different excitation properties. To obtain a computation that is robust towards such an excitation change, we have pursued a different idea than with the robust subspace residual from (4.7): The techniques developed in Section 1.3 are used to “normalize” Hankel matrices from different measurement setups under possibly different excitation properties. We can use a similar strategy here to make $\hat{\mathcal{H}}_{\text{ref}} = \hat{\mathcal{H}}_{\text{ref}}(\theta_0)$ and $\hat{\mathcal{H}}_{\text{test}} = \hat{\mathcal{H}}_{\text{test}}(\mathcal{Y}_N)$ compatible for a comparison under changing excitation by defining a residual with a normalized Hankel matrix from the test data

$$\tilde{\xi}(\theta_0, \mathcal{Y}_N) = \sqrt{N} \text{vec}(\hat{\mathcal{H}}_{\text{test}} \hat{\mathcal{Z}}_{\text{test}}^{\dagger} \hat{\mathcal{Z}}_{\text{ref}} - \hat{\mathcal{H}}_{\text{ref}}), \quad (4.9)$$

where the factors $\hat{\mathcal{Z}}_{\text{ref}}$ and $\hat{\mathcal{Z}}_{\text{test}}$ are defined from the truncated SVD at order $n$ of the juxtaposed Hankel matrix estimates

$$[\hat{\mathcal{H}}_{\text{ref}} \quad \hat{\mathcal{H}}_{\text{test}}] = [\hat{U}_1 \quad \hat{U}_2] \begin{bmatrix} \hat{D}_1 & 0 \\ 0 & \hat{D}_2 \end{bmatrix} \begin{bmatrix} \hat{V}_{1,\text{ref}}^T & \hat{V}_{1,\text{test}}^T \\ \hat{V}_{2,\text{ref}}^T & \hat{V}_{2,\text{test}}^T \end{bmatrix}, \quad (4.10)$$
as $\hat{Z}_{\text{ref}} = \hat{D}_1 \hat{V}_{1,\text{ref}}^T$, $\hat{Z}_{\text{test}} = \hat{D}_1 \hat{V}_{1,\text{test}}^T$. Besides the robustness of the residual towards excitation changes, we have proved that it fits into the asymptotic local approach framework, satisfying (4.2), and thus can be used for parametric change detection and isolation. Furthermore, the developed method can take into account the uncertainty inherent to both the reference and the tested data sets, unlike in previous approaches, which yields an improved noise robustness. The details are shown in [J23].

[J7] M. Döhler and L. Mevel, “Subspace-based fault detection robust to changes in the noise covariances,” *Automatica*, vol. 49, no. 9, 2734–2743, 2013

[J8] M. Döhler, L. Mevel, and F. Hille, “Subspace-based damage detection under changes in the ambient excitation statistics,” *Mechanical Systems and Signal Processing*, vol. 45, no. 1, 207–224, 2014

[J23] S. Greś, M. Döhler, P. Andersen, and L. Mevel, “Subspace-based Mahalanobis damage detection robust to changes in excitation covariance,” *Structural Control and Health Monitoring*, vol. 28, no. 8, e2760, 2021

[CA2] M. Döhler and L. Mevel, “Robust subspace based fault detection,” in *Proc. 18th IFAC World Congress*, Milan, Italy, 2011

[CA16] S. Greś, M. Döhler, and L. Mevel, “Hankel matrix-based Mahalanobis distance for fault detection robust towards changes in process noise covariance,” in *Proc. 19th IFAC Symposium on System Identification (SYSID)*, 2021

[CS33] M. Döhler and F. Hille, “Subspace-based damage detection on steel frame structure under changing excitation,” in *Proc. 32nd International Modal Analysis Conference (IMAC)*, Orlando, FL, USA, 2014

[CS64] S. Greś, M. Döhler, P. Andersen, L. Damkilde, and L. Mevel, “Hankel matrix normalization for robust damage detection,” in *Proc. 8th International Operational Modal Analysis Conference (IOMAC)*, Copenhagen, Denmark, 2019

### 4.2.3 Robustness to uncertainties in reference

In the definition of the residual (4.6) or (4.7), only the Hankel matrix estimate $\hat{\mathcal{H}}$ computed from the test data has been considered as a random variable in previous works [51]–[54], [J8], while the left null space $\mathcal{S}$ associated to the reference state is considered to be deterministic. Hence, both the residual distribution in (4.2) and the damage detection test (4.3) have been derived in previous works assuming that $\mathcal{S}$ is deterministic.

However, this is not a realistic assumption. In practice, in the absence of a reference model, only an estimate $\hat{\mathcal{S}}$ of $\mathcal{S}$ is usually available, which is computed from a Hankel matrix estimate $\hat{\mathcal{H}}_{\text{ref}}$ in the reference state. By nature, this Hankel matrix estimate is afflicted with statistical uncertainty and so is $\hat{\mathcal{S}}$. The related uncertainty may not be negligible, in particular for cases when the reference data set for the computation of $\hat{\mathcal{H}}_{\text{ref}}$ is relatively short.
Consequently, the residual mean in the reference state is not exactly zero anymore under small errors in $S$ and introduces bias, and the residual distribution in (4.2) may be incorrect. This can lead to an incorrect behavior of the damage detection test and to unwanted false alarms. To circumvent this problem, the uncertainty associated to the reference state is developed and included in the residual covariance in this line of work. In this way, the bias due to the error in identifying $S$ is removed by considering it as a random variable, so instead of having a bias the covariance is adapted appropriately.

Hence, instead of (4.6), the residual is considered as a function of two random variables

$$ \tilde{\zeta}(\theta_0, \mathcal{Y}_N) \stackrel{\text{def}}{=} \sqrt{N} \text{vec}(\hat{S}^T \hat{\mathcal{H}}), \quad (4.11) $$

where $N$ is the length of the test data to obtain $\hat{\mathcal{H}}$, and assume that $M$ is the length of the data in the reference state that is used to obtain $\hat{S}$. Then, the analysis of residual distribution yields its asymptotic covariance estimate

$$ \hat{\Sigma} = \frac{N}{M}\hat{\Sigma}_1 + \hat{\Sigma}_2, \quad (4.12) $$

where $\Sigma_1$ is related to the uncertainty in $\hat{S}$, and $\Sigma_2$ is related to the uncertainty of the current test data in $\hat{\mathcal{H}}$. The covariance contribution $\Sigma_2$ is exactly the same as the covariance of residual (4.6) when considering $S$ as deterministic, and $\frac{N}{M}\hat{\Sigma}_1$ is the additional covariance term due to uncertainties in $\hat{S}$, developed in detail in [J26]. Hence, when the data length $M$ to compute the reference matrix $\hat{S}$ is large with respect to the data length $N$ for computing $\hat{\mathcal{H}}$ from the current test data, the covariance contribution $\hat{\Sigma}_1$ is small, which is reflected in a small $\frac{N}{M}$ in (4.12). On the other side, when relatively few data $M$ are available in the reference state, the contribution of the uncertainty related to the reference matrix $\hat{S}$ is significant.

The correct consideration of uncertainties in the estimation of $S$ leads to the expected statistical properties of the associated test statistic (4.3), whose mean in the reference state is supposed to be the dimension of the parameterization $\theta$. In particular, this allows the definition of theoretical instead of empirical thresholds for the new test, whereas the bias in the conventional test computation has led to an unpredictable deviation. This can be

<figure><img src="image_19.png" /><figcaption>Figure 4.3 – Monte Carlo simulation on mass-spring chain with conventional ($\Sigma_1 = 0$) and new covariance computation, with $N = 100,000$. Left: Mean of test statistic in reference state for different $M$. Right: Probability of detection for small damages with $M = 50,000$.</figcaption></figure>
illustrated in Monte Carlo simulations on a mass-spring chain: In Fig. 4.3 (left), the mean of the test statistic is shown without considering the uncertainties in the reference (conventional computation) and with considering them (new computation) for different data lengths of the reference data. Indeed the test mean deviates significantly from the theoretical value when the uncertainties in the reference are not considered, whereas it matches very well when they are considered. Moreover, the correct consideration of the uncertainties has a significant impact on the test performance: In Fig. 4.3 (right), the probability of detection (POD) is shown for the same probability of false alarms in the conventional and new test computations. It is higher when correctly considering the uncertainties in the reference.

[J26] E. Viefhues, M. Döhler, F. Hille, and L. Mevel, “Statistical subspace-based damage detection with estimated reference,” *Mechanical Systems and Signal Processing*, vol. 164, 2021–2034, 2022

[CA11] E. Viefhues, M. Döhler, F. Hille, and L. Mevel, “Asymptotic analysis of subspace-based data-driven residual for fault detection with uncertain reference,” in *Proc. 10th IFAC Symposium on Fault Detection, Diagnosis and Safety of Technical Processes (SAFEPROCESS)*, Warsaw, Poland, 2018

[CS51] E. Viefhues, M. Döhler, F. Hille, and L. Mevel, “Stochastic subspace-based damage detection with uncertainty in the reference null space,” in *Proc. 11th International Workshop on Structural Health Monitoring (IWSHM)*, Stanford, CA, USA, 2017

### 4.2.4 Robustness to external nuisance

It is well-known that the dynamic system properties not only change with damage, but also with environmental nuisance, such as the temperature. Let $P$ be the external parameter that describes the nuisance. The described damage detection method works under the premise that the orthogonality between the left reference null space $S$ and the Hankel matrix estimate $\hat{H}$ in residual (4.6) (or (4.7), (4.11)) is only lost when damage appears, which means that it is only meaningful when the test data is obtained under the same external parameter as $S$. If the Hankel matrix is computed from test data under a different external parameter than the one of $S$, the residual mean is not zero anymore in the reference state, leading to false alarms of the damage detection test (4.3). In this case, a left reference null space would be required that is appropriate for the environmental conditions of the test data, which would avoid such false alarms and enables a test setup robust to environmental variations.

To formalize the damage detection problem in this situation, it is supposed that the test data is obtained under a known external parameter $P$ that is possibly different from a set of reference parameters $\{P_1, P_2, \dots, P_u\}$ at which measurements in the reference state of the structure are available. Based on these reference measurements, it is the goal to develop an appropriate damage detection test for data under parameter $P$, using model interpolation for linear parameter varying (LPV) systems [57].

Instead of one Hankel matrix $\hat{H}_{\text{ref}}$ in the reference state on which the left null space $\hat{S}$ is
computed, several reference Hankel matrices $\hat{\mathcal{H}}_j$ are available now under the different external parameters $P_j, j = 1, \dots, u$ that are in general different from $P$. Based on these references at values $\{P_1, P_2, \dots, P_u\}$, an appropriate left reference null space $S(P)$ under the current parameter $P$ of the testing state is created. We pursue two strategies: *Strategy 1* is based on output interpolation at the reference points [57], where a weighting function $\rho_j(P)$ is defined, e.g., a bell-shaped (Gaussian) function centered at $P$, and the weights $\rho_j(P), j = 1, \dots, u$ are the values of this function at $P_j$ after normalization such that $\sum \rho_j(P) = 1$. We have shown that this is equivalent to interpolating directly the Hankel matrices in the reference state, with

$$ \underline{\mathcal{H}}(P) = \sum_{j=1}^{u} \rho_j^2(P) \mathcal{H}_j, \qquad (4.13) $$

and the left reference null space $S(P)$ is computed on $\underline{\mathcal{H}}(P)$. *Strategy 2* is inspired by this interpolation approach, but instead of the Hankel matrices the modal parameters are interpolated. Let $\theta_j$ be the vector of the identified modal parameters at $P_j$, then

$$ \underline{\theta}(P) = \sum_{j=1}^{u} \rho_j(P) \theta_j. \qquad (4.14) $$

Based on $\underline{\theta}(P)$, the system matrices $\underline{A}(P)$ and $\underline{C}(P)$ and the respective observability matrix $\underline{Q}(P)$ can be built (see (1.5)), whose left null space $S(P)$ is obtained.

The obtained $S(P)$ is then confronted to current data $\mathcal{Y}_N(P)$ in residual (4.11). The covariance of this residual is thoroughly evaluated, accounting for the uncertainties related to the different reference datasets, and used in the detection tests.

Within the developed strategies, *strategy 1* is the simplest and the most straightforward to implement, as it is purely data-driven and modal analysis is not required. *Strategy 2* requires in addition modal parameter estimation, but only once from the reference measurements. Modal analysis in the testing state is never required. Using their interpolation to construct an appropriate left reference null space, strategy 2 considers more directly the well-known effect of environmental variation on the modal parameters than on the related Hankel matrices in strategy 1, and should be more robust to noise.

These strategies use the knowledge of the external parameter in addition to the vibration measurements (e.g., the temperature needs to be measured). This way, they have the potential to define a more precise reference for the current environmental condition than other approaches that mix data blindly from different but unknown environmental conditions for the definition of a reference. In particular, in the context of subspace-based damage detection under environmental variation, a previously developed approach is the *average approach* from [53], where the left reference null space is computed on the mean of the Hankel matrices $\frac{1}{u} \sum_{j=1}^{u} \mathcal{H}_j$ that are obtained from measurements in the reference state under different environmental conditions. While this approach provides some robustness in practice when the external parameter is unknown, it is not optimal particularly when the underlying external parameter is quite far from the mean.

This is illustrated on a mass-spring chain simulation, where a temperature dependent stiffness was defined. In Fig. 4.4, test values are shown for data in healthy and damaged
<figure><img src="image_20.png" /><figcaption>Figure 4.4 – Test values at different testing temperatures computed with the two developed strategies: Strategy 1 (S(H), top left), Strategy 2 (S(O), top right), and comparison to the average approach (bottom).</figcaption></figure>

states at different temperatures, where damage was a small stiffness reduction of 1.5% in one of the springs. Reference datasets were generated at temperatures {0, 5, 10, 15, 20}, and the threshold was defined based on these reference datasets. Strategy 1 leads to more false alarms and a less clearer distinction between healthy and damaged states than strategy 2, whereas the previous average approach fails to separate healthy and damaged states for the small damage. The developed methods were also applied to data of a test bridge in an outdoor lab of BAM, where Strategy 2 also led to the fewest false alarms even when only few reference datasets were used.

[CA13] E. Viefhues, M. Döhler, F. Hille, and L. Mevel, "Fault detection for linear parameter varying systems under changes in the process noise covariance," in *Proc. 21st IFAC World Congress*, 2020

[CS62] E. Viefhues, M. Döhler, Q. Zhang, F. Hille, and L. Mevel, "Subspace-based damage detection with rejection of the temperature effect and uncertainty in the reference," in *Proc. 8th International Operational Modal Analysis Conference (IOMAC)*, Copenhagen, Denmark, 2019
[CS71] E. Viefhues, M. Döhler, P. Simon, R. Herrmann, F. Hille, and L. Mevel, “Stochastic subspace-based damage detection of a temperature affected beam structure,” in *10th International Conference on Structural Health Monitoring of Intelligent Infrastructure (SHMII)*, 2021

### 4.2.5 Localization: robustness to over-parametrization

The robustness issues of the previous sections concern the definition and the statistical properties of the residual. This section focuses on its evaluation with respect to the parameterization $\theta$ with the particular purpose of damage (or in general *change*) localization. Damage is defined as a change in structural parameters, which are typically related to structural properties of a finite element model like material properties or cross-section values, prestressing forces, support conditions, stiffness of subcomponents, geometric properties and parameters that describe the system connectivity or local mass properties. The sensitivity matrix $\mathcal{J}$ in (4.2) is computed analytically with the help of the finite element model.

The main idea of damage localization is to test the hypotheses from Eq. (4.1) individually for each parameter component in $\theta$, and to localize damage by selecting the parameters that exhibit a test response beyond a prescribed threshold. In previous works [52], [54], the *direct test* (also called *sensitivity test*) in (4.5) has been used for this task. The direct test is straightforward in its implementation but neglects possible changes in the untested parameters, which can lead to significant test reactions for unchanged parameters beyond the prescribed threshold (false alarms) when parameter sensitivities are not orthogonal. Instead, we have proposed to use the minmax test [58] that avoids this problem by appropriately considering the untested parameters: When testing parameter $\theta_i$ for a change, define $\mathcal{J}_i$ as the corresponding column of $\mathcal{J}$ and denote the remaining columns by $\mathcal{J}_{\bar{\tau}}$. Define the partitioned Fisher information matrix

$$ \begin{bmatrix} F_{i,i} & F_{i,\bar{\tau}} \\ F_{\bar{\tau},i} & F_{\bar{\tau},\bar{\tau}} \end{bmatrix} = \begin{bmatrix} \mathcal{J}_i^T \Sigma^{-1} \mathcal{J}_i & \mathcal{J}_i^T \Sigma^{-1} \mathcal{J}_{\bar{\tau}} \\ \mathcal{J}_{\bar{\tau}}^T \Sigma^{-1} \mathcal{J}_i & \mathcal{J}_{\bar{\tau}}^T \Sigma^{-1} \mathcal{J}_{\bar{\tau}} \end{bmatrix}. \quad (4.15) $$

Then, define the partial residuals $\zeta_i \stackrel{\text{def}}{=} \mathcal{J}_i^T \Sigma^{-1} \zeta$ and $\zeta_{\bar{\tau}} \stackrel{\text{def}}{=} \mathcal{J}_{\bar{\tau}}^T \Sigma^{-1} \zeta$, the robust residual $\zeta_i^* \stackrel{\text{def}}{=} \zeta_i - F_{i,\bar{\tau}} F_{\bar{\tau},\bar{\tau}}^{-1} \zeta_{\bar{\tau}}$ and $F_i^* \stackrel{\text{def}}{=} F_{i,i} - F_{i,\bar{\tau}} F_{\bar{\tau},\bar{\tau}}^{-1} F_{\bar{\tau},i}$. With these projections, the mean of the robust residual $\zeta_i^*$ is only sensitive to changes $\delta_i$, but blind to $\delta_{\bar{\tau}}$, and the corresponding GLR

<figure><img src="image_21.png" /><figcaption>Figure 4.5 – Direct tests (left) and minmax tests (right) for six stiffness parameters of mass-spring chain with 5% damage in elements 4 and 6.</figcaption></figure>
test statistic for $\delta_i = 0$ against $\delta_i \neq 0$ turns out to be

$$t_{\text{mm}}^{i} = \zeta_{i}^{*T} F_{i}^{*-1} \zeta_{i}^{*}, \qquad (4.16)$$

denoted as the minmax test. As an example, Fig. 4.5 shows the test values for a mass-spring chain with six elements, where damage is introduced by a 5% stiffness decrease in elements 4 and 6. While the direct test shows also strong test reactions for elements 3 and 5, the minmax test reacts as expected.

In the damage localization tests, the parameter sensitivities need to be sufficiently different to be able to distinguish individual parameter changes. The sensitivity matrix $\mathcal{J}$ needs to be of full column rank, meaning changes in the damage-sensitive feature can be clearly related to the respective structural parameters. Due to the possible large number of structural parameters in FE models, this requirement is generally not given. Therefore, the parameters with similar sensitivity vectors have to be clustered prior to the damage localization. In the minmax test this is also a technical requirement for the projections that make the test blind to parameter changes in untested elements. Previously, k-means clustering was proposed [52], [54], which however depends strongly on a random initialization and does not necessarily cluster close elements satisfactorily. We propose a hierarchical clustering strategy, namely complete-linkage clustering [59], and show that it is well adapted to our damage localization approach.

With these developments that also include a revisited model-based sensitivity computation, the theoretical framework for statistical sensitivity-based damage localization was developed into a working damage localization method that is more robust for real applications. The resulting method offers a flexible and generic framework for damage localization, taking into account both physical model information and data-based uncertainties. The method was successfully applied to experimental data for the first time, namely from different damage cases of the Yellow Frame, a 3D steel frame, at the University of British Columbia, with some results shown in Fig. 4.6. To introduce damage, some of the 32 braces were removed,

<figure><img src="image_22.png" /><figcaption>Figure 4.6 – Damage localization tests for 32 brace elements. Left: damage scenario with removal of braces 21 and 23, right: damage scenario with removal of braces 2, 4, 18 and 20.</figcaption></figure>
and the tests were computed with respect to the stiffness parameterization of the brace elements. We could illustrate that all three ingredients (minmax tests, model-based sensitivity computation, hierarchical clustering) were necessary for a successful damage localization.

[J18] S. Allahdadian, M. Döhler, C. Ventura, and L. Mevel, “Towards robust statistical damage localization via model-based sensitivity clustering,” *Mechanical Systems and Signal Processing*, vol. 134, 106341, 2019

[CS55] S. Allahdadian, M. Döhler, C. Ventura, and L. Mevel, “Damage localization of a real structure using the statistical subspace damage localization method,” in *Proc. 11th International Workshop on Structural Health Monitoring*, Stanford, CA, USA, 2017

### 4.2.6 Numerical robustness of test computation

The damage detection and localization tests (4.3), (4.5) and (4.16) require a number of matrix inversions, which may be numerically critical due to possible ill-conditioning of the sensitivity and covariance matrices. The goal is to keep the number of matrix inversions to a minimum in the numerical computations of the test statistics, and to keep the sizes of the involved matrices small for an efficient computation.

The developments for numerically robust test computations are based on:

* The structure of covariance estimate $\hat{\Sigma}$ is exploited, similarly as in Section 2.2. In practice, the number of data blocks for its estimation is often lower than its dimension, leading to the low-rank decomposition $\hat{\Sigma} = \hat{\Sigma}^{1/2}(\hat{\Sigma}^{1/2})^T$ with an often much smaller matrix $\hat{\Sigma}^{1/2}$. This decomposition can (and should) also use truncation of small singular values if $\hat{\Sigma}$ is badly conditioned. In theory, $\hat{\Sigma}$ needs to be of full rank and its inverse is required in the tests. However, it can be shown that its truncated version still preserves the Gaussian properties of the residual, and the pseudoinverse can be used instead [J7]; hence only $(\hat{\Sigma}^{1/2})^\dagger$ is computed.

* Use of QR decompositions: In the global test, the thin QR decomposition $(\hat{\Sigma}^{1/2})^\dagger \hat{J} = QR$ is used, yielding an efficient computation $t_{\text{global}} = \alpha^T \alpha$ with $\alpha = Q^T(\hat{\Sigma}^{1/2})^\dagger \zeta$ in (4.3). In the minmax test for element $i$, the thin QR decomposition

$$ (\Sigma^{1/2})^{\dagger} [\mathcal{J}_{\bar{i}} \ \mathcal{J}_{i}] = [Q_{\bar{i}} \ Q_{i}] \begin{bmatrix} R_{\bar{i}\bar{i}} & R_{\bar{i}i} \\ 0 & R_{ii} \end{bmatrix} \quad (4.17) $$

is used, yielding $t_{\text{mm}}^i = \alpha^T \alpha$ with $\alpha = Q_i^T(\hat{\Sigma}^{1/2})^\dagger \zeta$ in (4.16). Furthermore, these QR decompositions can be computed iteratively at little cost when a large number of elements $i$ is tested.

It has turned out that these strategies are quite essential for a robust test computation on real world data. They ensure a numerically stable computation, and also speed up the computations significantly when the residual dimension is high.
[J13] M. Döhler, L. Mevel, and Q. Zhang, “Fault detection, isolation and quantification from Gaussian residuals with application to structural damage diagnosis,” *Annual Reviews in Control*, vol. 42, 244–256, 2016

[CA6] M. Döhler, L. Mevel, and F. Hille, “Efficient computation of minmax tests for fault isolation and their application to structural damage localization,” in *Proc. 19th IFAC World Congress*, Cape Town, South Africa, 2014

## 4.3 Kalman filter-based additive change detection and isolation

The asymptotic local approach framework (Section 4.2) relies on a CLT for the damage-sensitive feature in (4.2). An entire dataset of length $N$ is evaluated, and the used feature may be of considerable size, such as the subspace-based residual in (4.6). An alternative is provided by *additive change detection*, where the change in the system is modelled by additive terms in the state and output equations of the type

$$x_{k+1} = A_0 x_k + \Psi_k \tilde{\delta} + w_k \quad (4.18a)$$

$$y_k = C_0 x_k + \Phi_k \tilde{\delta} + v_k \quad (4.18b)$$

where $A_0$ and $C_0$ are the system matrices in the reference state corresponding to $\theta_0$, $\Psi_k$ and $\Phi_k$ are fault profile matrices that indicate how a parameter change affects the states and the outputs, and $\tilde{\delta} = \theta - \theta_0$. In particular, this modelling allows a sample-by-sample analysis without the need of a CLT, the handling of relatively small scale matrices of dimensions of only the number of outputs $r$ or of the system order $n$, and the application of classical techniques for additive change detection [55].

However, our parameter changes are not additive but multiplicative, i.e., they affect both the system matrices (eigenvalues and eigenvectors) and the states. To arrive at a formulation like (4.18), we propose a first-order perturbation analysis under the assumption of a small parameter change $\theta = \theta_0 + \varepsilon\theta_1$, $\tilde{\delta} = \varepsilon\theta_1$. Then, accordingly, for $A, C$ and $x_k$ under system parameter $\theta$,

$$A = A_0 + \varepsilon A_1, \quad (4.19a)$$

$$C = C_0 + \varepsilon C_1, \quad (4.19b)$$

$$x_k = x_k^0 + \varepsilon x_k^0, \quad (4.19c)$$

where $x_k^0 = \hat{x}_{k|k-1}$ is the one-step ahead state prediction from a Kalman filter based on the nominal system with $A_0$ and $C_0$. Plugging these approximations into the state-space model (1.2), and making explicit the link between changes in the eigenstructure of $(A, C)$ and changes in the parameterization $\theta$ through the respective sensitivities yields an approximation by model (4.18), where $\Psi_k$ and $\Phi_k$ are filled with known sequences (in particular based on the Kalman filter state predictions $x_k^0$).
Then, the mean of the innovations $\zeta_k = y_k - C_0 x_k^0$ can be directly related to parameter change $\tilde{\delta}$ [55], which is the most important result of these developments. Under the assumption of Gaussian noises, it holds

$$ \zeta_k \approx \begin{cases} \mathcal{N}(0, \Sigma_e) & \text{under } \mathbf{H}_0 \\ \mathcal{N}(\mathcal{J}_k \tilde{\delta}, \Sigma_e) & \text{under } \mathbf{H}_1 \end{cases} \qquad (4.20) $$

where $\mathcal{J}_k = C_0\Gamma_k + \Phi_k$ and $\Gamma_k$ is recursively computed as

$$ \Gamma_{k+1} = C_0(I_n - K C_0)\Gamma_k + \Psi_k - A_0 K \Phi_k, $$

starting with $\Gamma_0 = 0$; $K$ is the Kalman gain obtained from the nominal model, and $\Sigma_e$ is the innovation covariance. Thanks to (4.20) we have returned to the Gaussian linear framework, where the GLR tests for detection or localization can be applied analogously as introduced in Sections 4.2.1 and 4.2.5. For example, the global test writes

$$ s_{\text{global}} = \beta^T \Omega^{-1} \beta, \quad \text{where} \quad \beta = \sum_{k=1}^{N} \mathcal{J}_k^T \Sigma_e^{-1} \zeta_k, \quad \Omega = \sum_{k=1}^{N} \mathcal{J}_k^T \Sigma_e^{-1} \mathcal{J}_k \qquad (4.21) $$

In this formulation, the involved matrices are in general smaller and better conditioned than with the subspace-based residual in previous sections. In particular, the innovation covariance matrix is only of size of the number of sensors. It should be easily invertible without the issues described in Section 4.2.6. Moreover, a sample-wise evaluation is possible. On the other side, the noise covariance matrices of the system need to be estimated to set up the Kalman filter, which is not necessary when evaluating the subspace-based residual. The Kalman-based framework has been successfully tested on experimental data for damage detection and on simulated data for damage localization.

[CA8] M. Döhler, Q. Zhang, and L. Mevel, “Vibration monitoring by eigenstructure change detection based on perturbation analysis,” in *Proc. 17th IFAC Symposium on System Identification (SYSID)*, Beijing, China, 2015

[CA9] M. Döhler, Q. Zhang, and L. Mevel, “Change detection and isolation in mechanical system parameters based on perturbation analysis,” in *Proc. 20th IFAC World Congress*, Toulouse, France, 2017

[CA15] M. Döhler, Q. Zhang, and L. Mevel, “Damage localization in mechanical systems by Lasso regression,” in *Proc. 19th IFAC Symposium on System Identification (SYSID)*, 2021

[CA17] M. Döhler, Q. Zhang, and L. Mevel, “Kalman predictor subspace residual for mechanical system damage detection,” in *Proc. 11th IFAC Symposium on Fault Detection, Diagnosis and Safety of Technical Processes (SAFEPROCESS)*, Pafos, Cyprus, 2022
## 4.4 Conclusions and outlook

The contributions of this chapter have two different levels of maturity. In Section 4.2, theory is developed to make damage diagnosis methods in the asymptotic local approach framework robust under typical application conditions. This results in statistical data-driven algorithms that are automated, fast and numerically robust, and thus suitable for online monitoring. They are robust under naturally changing excitation properties and environmental nuisance (noise, temperature changes). The strongly improved applicability of these methods on real data also led to industrial transfer. The robust implementation of the damage detection method (Sections 4.2.2 and 4.2.6) has been transferred to ARTEMIS Modal. Furthermore, the first application of damage localization became possible on experimental data (beyond numerical simulations) with this method class. Since damage localization (and quantification) requires in addition a physical model, it is a more complex problem than detection, and the transition still needs to be made from laboratory experiments to real-world case studies with complex and large-scale structures under realistic damage. This certainly requires addressing some further numerical and statistical challenges, related e.g. to optimal parameter clustering, or to consider model uncertainties.

Section 4.3 is more exploratory, going back to the foundations of the statistical damage diagnosis frameworks. There we have developed an alternative Kalman-filter based approach in a more theoretical contribution, whose full potential for monitoring still needs to be explored, in particular for real-time monitoring, and more general for physical model-based data analysis for digital twins.
# Damage diagnosis performance and reliability

## 5.1 Introduction

With a statistical method for damage diagnosis at hand, the question arises to which kinds of structural damages can actually be detected or localized, or, more precisely, how big must a damage (a parameter change) be so that it can be reliably diagnosed with the method under the typical data-based uncertainties. This way, the performance of the considered damage diagnosis method can be assessed even before damage occurs, which is an essential information for stakeholders or decision makers to deploy SHM in practice. Surprisingly, this is a question that has hardly been considered in the SHM context in the literature [60], and appropriate tools like the probability of detection (POD) curves are typically used in the NDT (non-destructive testing) context where they are obtained based on data from damaged specimen, but hardly in SHM where usually no data from the damaged state are available.

An important feature of the statistical frameworks for damage detection and localization from Chapter 4 is the characterization of the probability distributions of the damage detection and localization tests for both reference and damaged states. Based on this framework, we are able to define and to analyze the damage detectability and localizability for the respective methods. No data from the damaged states are required, but information from the physical model (in the reference state only) is used instead via the model-based parameter sensitivities. With the POD of a damage diagnosis method, a step forward can be made to bridge the gap between SHM-based damage diagnosis and the actual structural reliability and performance assessment: the information on the health state (from data) and the POD of expected damages during the service life of a structure can be used to update the system reliability based on a deterioration model. In combination with expected costs, the value of information of the data can also be evaluated in the design process of an SHM system.
## Contributions

My contributions are methods for the evaluation of diagnosis performance and structural reliability within the local approach framework:

* Definition and evaluation of the damage detectability and localizability. This allows to assess the performance of damage detection and localization methods in this framework, and yields indicators that allow the optimization of the performance, e.g., with optimal sensor placement (Section 5.2).

Context: PhD of Alexander Mendler at University of British Columbia, Canada (supervision together with Carlos Ventura).

* Bayesian update of structural system reliability and evaluation of the value of information in a further step towards performance assessment in this framework, by exploiting the detectability together with deterioration and cost models, respectively (Section 5.3).

Context: Collaboration with Sebastian Thöns and PhD student Lijia Long at BAM.

## 5.2 Detectability and localizability

The actual performance of a damage detection or localization method is hardly analyzed in general, and it is often unclear which damages can actually be detected or localized with a given detection/localization method for the given measurement setup. For this, the notions of detectability and localizability have to be established first. To this end, damage is defined as a change in system parameter vector $\theta$, as in Chapter 4. For simplicity of notation, we suppose individual parameter changes.

**Detectability.** We define a damage in parameter component $\theta_i$ as detectable with methods defined by (4.1)–(4.3) when the corresponding change from the reference value $\theta_i^0 - \theta_i$ yields a non-centrality parameter $\lambda$ of the test that exceeds a minimum value $\lambda \ge \lambda_{min}$. This means that the distributions of the test statistic in the reference and damaged states (as depicted in Fig. 4.1) need to be sufficiently separated. The value $\lambda_{min}$ relates to reliability considerations, and can be uniquely defined with the $\chi^2$ distribution properties based on the acceptable probability of false alarms (PFA) and the targeted probability of detection (POD) for which a damage is considered as detectable. Once $\lambda_{min}$ is set, the required minimal parameter change $\theta_i^0 - \theta_i$ that yields $\lambda \ge \lambda_{min}$ can be computed from property (4.4) of the non-centrality parameter. Plugging in $\delta = \sqrt{N}(\theta_i^0 - \theta_i)$ results in

$$ \theta_i^0 - \theta_i = \sqrt{\frac{\lambda_{min}}{N \cdot F_{ii}}}, \qquad (5.1) $$

where $F_{ii} = \mathcal{J}_i^T \Sigma^{-1} \mathcal{J}_i$ the i-th entry on the diagonal of the Fisher information $F$. Hence, the developments for the minimum detectable damage consider the statistical variability of the damage-sensitive feature, its link to structural design parameters through FE models (both contained in $F$), and the desired reliability of the damage detection results. The reliability considerations are incorporated including a user-defined PFA and POD, which
<figure><img src="image_23.png" /><figcaption>Figure 5.1 – Empirical $\chi^2$ distribution from Monte Carlo experiment with experimental data (training and testing state) while monitoring mass changes. Results for an 5% extra mass on Segment 5 (left) and on Segment 8 (right) of a beam.</figcaption></figure>

can, for example, be chosen based on design standards. Hence, the following aspects directly affect the minimum detectable damage:

(a) Damage-sensitive feature and its sensitivity

(b) Signal-to-noise ratio

(c) Measurement duration

(d) Reliability requirements for the damage detection results (i.e., PFA in the reference state, and the desired POD)

Maximum detectability is given for a high sensitivity of the feature, a high signal-to-noise ratio, a long measurement duration, and low requirements regarding the reliability, as can be directly inferred from (5.1).

The findings were validated on Monte Carlo simulations as well as on experimental data of a beam in the lab, where a change was introduced by adding masses. Using several datasets in the reference state and under a 5% mass change in selected sections of the beam, we could produce histograms of the test statistic distribution in both cases. Since the added mass was fixed, the measurement duration *N* was adapted in (5.1) to yield a $\lambda_{min}$ corresponding to 95% POD for the 5% mass change. An example is given in Fig. 5.1, showing that the empirically obtained PODs correspond very well.

**Localizability.** Damage localization is more challenging than the detection, first of all because not only one global test is computed but one test for each parameter. Due to a possibly large number of parameters in the model clustering may be required (see Section 4.2.5). If the clustering is not “perfect”, also undamaged elements might react in the tests, causing false alarms. Furthermore, the number of clusters is a natural limit of the localization capabilities, defining the possible resolution of the localization, and is therefore also an element to be considered when defining localizability. Hence, the definition of damage localizability is not straightforward, and we provide a framework based on the minimum detectable parameter changes of the tested elements, possible false alarms in undamaged elements, as well as the achievable damage localization resolution. The minimum localizable parameter change is derived analogously as in (5.1) when replacing $F_{ii}$ by $F_i^*$ (see Section 4.2.5). The analysis of the false alarm rate is quite original, as it can be derived from the Fisher information after clustering. While the mean of the robust residual $\zeta_i^*$ is in theory zero for undamaged elements
thanks to the involved projections, this is not ensured when plugging in the corresponding Fisher information matrices after the clustering process. From this quite technical analysis (see [CA12] for details) the residual mean and the resulting test reaction can be evaluated for undamaged elements, so false alarms can be known beforehand. Finally, all three *localizability criteria* (detectable parameter changes, false alarms, resolution) can be combined into a performance indicator for evaluation. As a side-product of this performance analysis, the localizability can be evaluated for different numbers of clusters (or different cluster configurations in general), in order to choose the configuration that optimizes the performance of the damage localization.

**Optimal sensor placement.** Different sensor configurations will influence the damage detectability and localizability analyzed before. Whereas the objective of traditional sensor placement methods is to maximize the identifiability and distinguishability of modes, or to reduce estimation uncertainty, they hardly target method-specific performance criteria that are tailored to damage diagnosis [61]. None of the existing criteria appear to consider the detectability of changes in selected structural components, which is our focus for optimal sensor placement. Considering that optimal damage detectability in all components with a limited number of sensors is unrealistic, the goal is to develop an optimization criterion that targets a prioritized set of damage scenarios. Such damage hotspots are typically known from the operational evaluation of the structure and should be considered during sensor placement. This way, the damage detectability is evaluated for each of the chosen parameters, and the sensor placement is chosen that optimizes the detectability for all chosen parameters. When the number of possible sensor placements is large, we have proposed the use of an adapted genetic algorithm. The result is a sensor layout that is more sensitive to anticipated

<figure><img src="image_24.png" /><figcaption>Figure 5.2 – Optimizing the test sensor layout using for damage (added mass) in element 8, and comparing the test statistics for real data.</figcaption></figure>
damage scenarios. Based on a laboratory steel beam, it is demonstrated that the findings
based on simulated vibration records can be transferred to real structures, despite the noisy
measurement environment and possibly different excitation levels. In Fig. 5.2, the histograms
of the test in reference and damaged state are shown based on experimental data for the sensor
configurations obtained from the developed method with the best (left), medium (middle)
and worst (right) performance. The performance is well reflected in the histograms.

[J24] A. Mendler, M. Döhler, and C. E. Ventura, “A reliability-based approach to determine the minimum detectable damage for statistical damage detection,” *Mechanical Systems and Signal Processing*, vol. 154, 107561, 2021

[J28] A. Mendler, M. Döhler, and C. E. Ventura, “Sensor placement with optimal damage detectability for statistical damage detection,” *Mechanical Systems and Signal Processing*, vol. 170, 108767, 2022

[J29] A. Mendler, M. Döhler, C. E. Ventura, and L. Mevel, “Localizability of damage with statistical tests and sensitivity-based parameter clusters,” *Mechanical Systems and Signal Processing*, vol. 204, 110783, 2023

[CA12] A. Mendler, M. Döhler, C. E. Ventura, and L. Mevel, “Clustering of redundant parameters for fault isolation with Gaussian residuals,” in *Proc. 21st IFAC World Congress*, 2020

[CS65] A. Mendler, S. Allahdadian, M. Döhler, L. Mevel, and C. Ventura, “Minimum detectable damage for stochastic subspace-based methods,” in *Proc. 8th International Operational Modal Analysis Conference (IOMAC)*, Copenhagen, Denmark, 2019

[CS68] A. Mendler, M. Döhler, C. E. Ventura, and L. Mevel, “The minimum detectable damage as an optimization criterion for performance-based sensor placement,” in *Proc. 38th International Modal Analysis Conference (IMAC)*, Houston, TX, USA, 2020

[CS72] A. Mendler, A. Cadoret, C. Freyssinet, M. Döhler, Y. Lecieux, L. Mevel, and C. Ventura, “Minimum localizable damage for stochastic subspace-based damage diagnosis,” in *10th International Conference on Structural Health Monitoring of Intelligent Infrastructure (SHMII)*, 2021

**5.3 Structural reliability updating and value of information**

The developments in the previous sections have shown that the performance of the damage
diagnosis can be evaluated based on data from the reference state (which is required to
estimate the covariance Σ) and the sensitivities J, under the small change hypothesis. In
particular, this allows to compute the probabilities of detection (POD) for any parameter
change based on (4.4) in a very simple way, without the need of estimating data or evaluating
a model in the damaged state. Results for the considered parameter change are accurate when errors introduced through the linearization in (4.2) are low. Such information is available thanks to the complete statistical characterization of the damage sensitive feature in both reference and damaged states, and is a particular strength of the local approach framework.

Probabilistic structural system performance models predict the deterioration of structural elements during the service life of a structure. Amongst others, they are used to evaluate the reliability of structural components or the entire structural system, e.g., via failure probabilities. Another application is in value of information (VOI) analyses, where they are used to predict future damages, and the difference between the monitoring cost and the saved cost by detecting those damages with an SHM system (the expected utility gain) corresponds to the value of SHM information. In both cases, the analyses are made using integrals over the space of damage scenarios that are described by the probabilistic performance model. When the reliability or VOI analysis are made in connection with the employed damage diagnosis method, then the POD needs to be evaluated for each realization of a damage scenario for the numerical evaluation of the integrals. For this task, the described framework for damage diagnosis is perfectly suited, since the POD can be easily computed with hardly any additional computational effort.

**Structural reliability updating.** It is our goal to incorporate information from monitoring when evaluating the reliability. Suppose that the system failure probability $P(F_s)$ is evaluated at a certain point in the service life of the structure where the system parameter $\theta$ is now a realization of random variable $\Theta$ that describes deterioration of the structure. When at this point the result of damage detection based on actual data becomes available, the failure probability $P(F_s)$ computed from the structural performance model can be updated with the damage detection result, since the POD for realizations of $\Theta$ are easily evaluated. Here we limit ourselves to the case where no damage is detected, denoted by the event “no indication”, $\bar{I}_s$. The respective conditional probability writes

$$P(F_s|\bar{I}_s) = \frac{P(F_s \cap \bar{I}_s)}{P(\bar{I}_s)} \qquad (5.2)$$

The probabilities in the numerator and denominator are obtained by numerical integration. E.g.,

$$P(\bar{I}_s) = \int_{\Omega_{\bar{I}_s}} (1 - P(I_s|\theta))p_{\Theta}(\theta)d\theta, \qquad (5.3)$$

where $p_{\Theta}$ is the probability density function of $\Theta$, $\Omega_{\bar{I}_s}$ is the integration domain defined by an appropriate limit state function, and $P(I_s|\theta)$ is the POD for a realization $\theta$.

**Value of information analysis.** The performance of the damage diagnosis method is an important part in the VOI analysis of an SHM system, as a higher performance will obviously lead to a higher VOI. Similarly as before, the structural performance model yields some $\Theta$ for a point in the service life. The POD related to realizations of $\Theta$ condition the actions (like “repair” or “no repair”) in a decision tree, updating the corresponding failure probabilities of the model and evaluating the utility gain of the monitoring. In a simple
simulation of a truss bridge girder we have evaluated the influence of different properties of the measurement system (number and position of sensors, measurement noise level, PFA) on the VOI to determine the configuration with the highest utility. For example, additional sensors or higher quality sensors with less noise increase cost, while they also increase the POD, so it depends on the cost modelling and the decision tree if they will increase the VOI or not. An example of such an analysis using the subspace-based residual (4.6) for damage detection is given in Fig. 5.3, where a relative VOI is evaluated for different numbers of sensors (left), sensor positions (middle) and measurement noise (right), when the SHM system is installed in a certain year during the service life.

<figure><img src="image_25.png" /><figcaption>Figure 5.3 – Value of information (y axis) of a truss with 29 elements for implementation of SHM system in different service years (x axis).</figcaption></figure>

[J16] S. Thöns, M. Döhler, and L. Long, “On damage detection system information for structural systems,” *Structural Engineering International*, vol. 28, no. 3, 255–268, 2018

[J27] L. Long, M. Döhler, and S. Thöns, “Determination of structural and damage detection system influencing parameters on the value of information,” *Structural Health Monitoring*, vol. 21, no. 1, 19–36, 2022

[CS23] S. Thöns and M. Döhler, “Structural reliability updating with stochastic subspace damage detection information,” in *Proc. 5th European Conference on Structural Control (EACS)*, Genoa, Italy, 2012

[CS41] M. Döhler and S. Thöns, “Efficient structural system reliability updating with subspace-based damage detection information,” in *Proc. 8th European Workshop on Structural Health Monitoring (EWSHM)*, Bilbao, Spain, 2016

[CS53] L. Long, S. Thöns, and M. Döhler, “Damage detection and deteriorating structural systems,” in *Proc. 11th International Workshop on Structural Health Monitoring (IWSHM)*, Stanford, CA, USA, 2017
[CS59] L. Long, S. Thöns, and M. Döhler, “The effects of deterioration models on the value of damage detection information,” in *6th International Symposium on Life-Cycle Civil Engineering (IALCCE)*, Ghent, Belgium, 2018

[CS56] L. Long, S. Thöns, and M. Döhler, “The effects of SHM system parameters on the value of damage detection information,” in *Proc. 9th European Workshop on Structural Health Monitoring (EWSHM)*, Manchester, UK, 2018

## 5.4 Conclusions and outlook

In this chapter, we have made a step forward to bridge the gap between SHM-based damage diagnosis and the actual structural reliability and performance assessment. Based on a combined data- and model-based characterization of the distribution properties of the damage diagnosis tests, a reliability-based analysis of the damage diagnosis performance has become possible. After being fed into structural deterioration models, this has led to methods for data-based updating of the structural reliability, as well as an analysis of the value of monitoring information.

A common point of all methods in Chapter 4 and this chapter is the use of sensitivities with respect to (structural) parameters for evaluation, under a small change hypothesis. However, their actual purpose is quite different: In Chapter 4 they are used to evaluate the damage-sensitive feature with respect to the parameters (in particular to do localization), whereas they are used to analyze performance of the damage diagnosis in a predictive sense in this chapter:

* What will be the performance of the method for some parameter change? or: How big does the parameter change have to be so that we will be able to detect it reliably? (Section 5.2)

* What will be the performance of the damage diagnosis during the service life for parameter changes according to a probabilistic deterioration model? How does the consideration of monitoring data then impact the failure probabilities during the service life, and will costs caused by repairs or structural failure actually be reduced when using monitoring-based maintenance? (Section 5.3)

The work on SHM performance assessment and optimization is still in its infancy. Although essential for many industry applications, almost no analytical methods are available for the rigorous performance analysis of damage diagnosis procedures before damage occurs, both in SHM and also in Non-Destructive Testing (NDT). Together with Alexander Mendler at TU München, current activities are related to establishing “predictive” POD curves for SHM and NDT [J31], and to develop the analysis for higher level damage diagnosis (localization, quantification).
Perspectives

The contributions of my work have led to methods with different levels of maturity for ap-
plications. The developed methods for system identification and the related uncertainty
quantification can be considered as quite mature for applications in vibration analysis at the
current stage. My contributions to statistical damage diagnosis have increased the robust-
ness of the methods under realistic application conditions. Therein, damage detection has
reached a relatively high level of maturity, whereas automated methods for higher level dam-
age diagnosis are the least developed so far. Furthermore, the methods for SHM performance
assessment and SHM-based reliability assessment are rather at the start of development.

Amongst others, remaining challenges lie in the automation and robustness of methods
to improve their industrial application level. Methods for higher level damage diagnosis and
performance analysis that work in connection with FE model-based information still have to
make their proof of concept for SHM on real-world data, involving complex structures and
models. The methodological numerical and statistical challenges related to such complex
applications need to be solved. However, real test cases of monitored structures before and
after actual damage events (beyond artificially introduced damage) are hardly available today.
A related challenge is the performance prediction of damage diagnosis methods, which is a
topic that is hardly considered currently, whereas the (proven) effectiveness of a damage
diagnosis method could help to convince decision-makers of the value of implementing a
monitoring system. Besides, being able to optimize its performance may improve the success
of future monitoring systems. Such developments for performance analysis, in particular for
damage localization and quantification, would also be beneficial to the development of digital
twins. Furthermore, going also into this direction is the problem of actual structural reliability
assessment and prediction based on monitoring data. Finally, new sensing technologies like
camera-based measurements will enable whole new monitoring opportunities, and at the
same time they will challenge the existing methods from modal analysis to damage diagnosis,
requiring adaptation and conception of new methods.

My future research project follows the previous activity and addresses some of these
challenges. Current and future work revolves around the following topics:

Uncertainty quantification and automated modal analysis

Methods for subspace-based system identification and uncertainty quantification (Chapters 1 and 2) are rather mature and have been subject of several technology transfers to industrial
partners, where the focus shifts from research to development.

A subject that remains on the theoretical side is related to the model order choice. While subspace-based system identification is consistent when the model order is chosen too high, the situation is different for the uncertainty quantification. With Szymon Greś, we have shown in very recent work that one of the first steps in the uncertainty propagation, i.e., the sensitivity computation of the SVD of the output covariance Hankel matrix is sensitive to the actual model order choice [CS79], [CS82]. This is particularly relevant for the analysis of stabilization diagrams, where the system is identified at many successive model orders. The subject requires further investigation, e.g., on how to optimize the model order choice when computing the SVD sensitivities, or to investigate the impact of an over-determined order on the actual accuracy of the results.

Concerning development, I am supervising research engineer Johann Priou in an ANR France Relance development project (2022–2024) for automated modal analysis and tracking. The goal is the development of a robust algorithm for automated extraction of the global modal parameters from the stabilization diagram, using the uncertainty information in a tailor-made clustering approach. Further objectives are the development of a robust tracking algorithm and environmental nuisance removal. First results are published in [CS74].

## Towards higher level damage diagnosis with complex structures

Methods for damage detection are data-driven and have reached some maturity for applications, as e.g. shown for the subspace-based damage detection on a progressive damage test on the S101 Bridge [J9]. Damage localization and quantification require in addition a physical model in order to link the changes in the data to physical changes of the structure for damage assessment, and are thus more complex. Automated methods for localization and quantification from Chapters 3 and 4, as well as the model-informed methods for damage detectability and localizability analysis in Chapter 5 have so far been validated on numerical simulations or applied on small-scale laboratory experiments like beams or frames, but not on real structures with complex models that are usually of large scale. This step is so far missing to ensure practical relevance and impact, and may require in addition the consideration of further statistical and numerical aspects, e.g., for an optimal clustering of the large number of model parameters (or a reduced-order modeling) with respect to the limited available information from sensor data, and consideration of both model-based and data-based uncertainties. Furthermore, questions arise towards the evaluation of the model accuracy that is actually necessary so that (sensitivity-based) damage localization and quantification approaches can work.

In this context, we are currently working together with Yann Lecieux at Nantes University on a laboratory case study on a mock-up of Saint-Nazaire Bridge. This is a cable-stayed bridge, where damage is introduced by cable failure. First results on damage localization with the DLV method from Section 3.2 are promising [CS83], and the goal is a full benchmark study.
## Vision-based vibration monitoring

New sensor technologies offer (near) full field measurements, which is a change of paradigm compared to measurements by a limited number of traditional sensors like accelerometers that are physically attached to the monitored structure. For example, the recent development of vision-based measurement technology through high speed cameras with associated motion extraction algorithms allows high resolution measurements, where essentially each pixel can serve as a sensor. These technologies may open new perspectives especially for damage diagnosis, since it offers much more spatial information. However, with the extremely high-dimensional data traditional techniques for system identification and damage diagnosis break down, which requires the conception of adequate methods not only from a numerical standpoint, but also from a statistical standpoint considering the highly increased spatial correlation of such data.

Currently, I am supervising PhD student Zhilei Luo (2021–2024) together with Qinghua Zhang at Inria and Vincent Baltazart at Université Gustave Eiffel, on the topic of robustness for vision-based vibration monitoring. In the thesis project, challenges related to typical application conditions in outdoor environments are addressed, like varying weather and light conditions that pose a problem for motion extraction methods. Furthermore, system identification with a huge number of sensors is addressed, based on PCA-based dimension reduction while keeping the high spatial resolution of the mode shapes. First results have been published in [CA18], [CS81].

## Kalman filter-based diagnosis and disturbance rejection

The Kalman filter-based framework in Section 4.3 provides an interesting alternative for damage detection, localization and possibly also quantification. So far, the developments are rather theoretical, and the full potential of the approach for damage diagnosis still needs to be explored, especially for real-time monitoring. A particular advantage of this framework for damage diagnosis is that – being built upon Kalman filter properties – features from filtering can be relatively easily included, like the ability to reject disturbances in a general setting [62]. Furthermore, the subject also connects to identification when it comes to the estimation of the noise covariance matrices from data that are required for the filters. A further important aspect to be explored is the suitability for digital twins, since the Kalman filter-based diagnosis framework offers a direct connection between the physical model of a structure and measurements in real time. This also relates to the more general question on how to combine data with physics-based models in the creation of digital twins.

## Damage diagnosis performance and reliability assessment

The developments of Chapter 5 provide a theoretical framework to predict damage detection performance under damage. The full potential of the approach still needs to be explored and validated on more complex structures. One goal is to predict full POD curves and their uncertainties accurately, where also model uncertainties should be included. Such a prediction can be used for the performance-based design of an SHM system, which should be explored in more detail. Similar results for higher level damage diagnosis (localization/quantification) are
currently missing [60], requiring first an accurate description of the probability distributions
of the damage localization/quantification results in dependence of actual changes in structural
model parameters. Such information should be relevant for SHM-based reliability assessment,
by providing an appropriate computational framework for the link with realistic deterioration
models, for example in a Bayesian updating context.

The overall goal of future research activities is the mathematical design of automated and
robust methods for a full SHM-based diagnosis of structures, exploiting data and physical
models in an optimal way and involving new monitoring technologies, on the way to future
intelligent infrastructure.
Publications

Articles in peer-reviewed international journals

[J1] M. Döhler, S. Kunis, and D. Potts, "Nonequispaced hyperbolic cross fast Fourier transform," *SIAM Journal on Numerical Analysis*, vol. 47, no. 6, 4415-4428, 2010.

[J2] M. Döhler and L. Mevel, "Fast multi-order computation of system matrices in subspace-based system identification," *Control Engineering Practice*, vol. 20, no. 9, 882-894, 2012.

[J3] M. Döhler and L. Mevel, "Modular subspace-based system identification from multi-setup measurements," *IEEE Transactions on Automatic Control*, vol. 57, no. 11, 2951-2956, 2012.

[J4] M. Döhler, X.-B. Lam, and L. Mevel, "Uncertainty quantification for modal parameters from stochastic subspace identification on multi-setup measurements," *Mechanical Systems and Signal Processing*, vol. 36, no. 2, 562-581, 2013.

[J5] M. Döhler and L. Mevel, "Efficient multi-order uncertainty computation for stochastic subspace identification," *Mechanical Systems and Signal Processing*, vol. 38, no. 2, 346-366, 2013.

[J6] M. Döhler, L. Marin, D. Bernal, and L. Mevel, "Statistical decision making for damage localization with stochastic load vectors," *Mechanical Systems and Signal Processing*, vol. 39, no. 1-2, 426-440, 2013.

[J7] M. Döhler and L. Mevel, "Subspace-based fault detection robust to changes in the noise covariances," *Automatica*, vol. 49, no. 9, 2734-2743, 2013.

[J8] M. Döhler, L. Mevel, and F. Hille, "Subspace-based damage detection under changes in the ambient excitation statistics," *Mechanical Systems and Signal Processing*, vol. 45, no. 1, 207-224, 2014.

[J9] M. Döhler, F. Hille, L. Mevel, and W. Rücker, "Structural health monitoring with statistical methods during progressive damage test of S101 Bridge," *Engineering Structures*, vol. 69, 183-193, 2014.

[J10] D. Bernal, M. Döhler, S. Mozaffari Kojidi, K. Kwan, and Y. Liu, "First mode damping ratios for buildings," *Earthquake Spectra*, vol. 31, no. 1, 367-381, 2015.
[J11] L. Marin, M. Döhler, D. Bernal, and L. Mevel, "Robust statistical damage localization with stochastic load vectors," *Structural Control and Health Monitoring*, vol. 22, no. 3, 557-573, 2015.

[J12] P. Mellinger, M. Döhler, and L. Mevel, "Variance estimation of modal parameters from output-only and input/output subspace-based system identification," *Journal of Sound and Vibration*, vol. 379, 1-27, 2016.

[J13] M. Döhler, L. Mevel, and Q. Zhang, "Fault detection, isolation and quantification from Gaussian residuals with application to structural damage diagnosis," *Annual Reviews in Control*, vol. 42, 244-256, 2016.

[J14] G. Gautier, L. Mevel, J.-M. Mencik, R. Serra, and M. Döhler, "Variance analysis for model updating with a finite element based subspace fitting approach," *Mechanical Systems and Signal Processing*, vol. 91, 142-156, 2017.

[J15] M. Bhuyan, M. Döhler, Y. Lecieux, L. Mevel, and F. Schoefs, "Statistical damage localization with stochastic load vectors using multiple mode sets," *Structural Health Monitoring*, vol. 16, no. 5, 518-535, 2017.

[J16] S. Thöns, M. Döhler, and L. Long, "On damage detection system information for structural systems," *Structural Engineering International*, vol. 28, no. 3, 255-268, 2018.

[J17] M. Bhuyan, G. Gautier, N. Le Touz, M. Döhler, F. Hille, J. Dumoulin, and L. Mevel, "Vibration-based damage localization with load vectors under temperature changes," *Structural Control and Health Monitoring*, vol. 26, no. 11, e2439, 2019.

[J18] S. Allahdadian, M. Döhler, C. Ventura, and L. Mevel, "Towards robust statistical damage localization via model-based sensitivity clustering," *Mechanical Systems and Signal Processing*, vol. 134, 106341, 2019.

[J19] S. Greś, M. Döhler, P. Andersen, and L. Mevel, "Kalman filter-based subspace identification for operational modal analysis under unmeasured periodic excitation," *Mechanical Systems and Signal Processing*, vol. 146, 106996, 2021.

[J20] S. Greś, M. Döhler, P. Andersen, and L. Mevel, "Uncertainty quantification for the modal phase collinearity of complex mode shapes," *Mechanical Systems and Signal Processing*, vol. 152, 107436, 2021.

[J21] S. Greś, M. Döhler, and L. Mevel, "Uncertainty quantification of the modal assurance criterion in operational modal analysis," *Mechanical Systems and Signal Processing*, vol. 152, 107457, 2021.

[J22] S. Greś, M. Döhler, and L. Mevel, "Statistical model-based optimization for damage extent quantification," *Mechanical Systems and Signal Processing*, vol. 160, 107894, 2021.

[J23] S. Greś, M. Döhler, P. Andersen, and L. Mevel, "Subspace-based Mahalanobis damage detection robust to changes in excitation covariance," *Structural Control and Health Monitoring*, vol. 28, no. 8, e2760, 2021.
[J24] A. Mendler, M. Döhler, and C. E. Ventura, "A reliability-based approach to determine the minimum detectable damage for statistical damage detection," *Mechanical Systems and Signal Processing*, vol. 154, 2021.

[J25] S. Greś, M. Döhler, N.-J. Jacobsen, and L. Mevel, "Uncertainty quantification of input matrices and transfer function in input/output subspace system identification," *Mechanical Systems and Signal Processing*, vol. 167, 2022.

[J26] E. Viefhues, M. Döhler, F. Hille, and L. Mevel, "Statistical subspace-based damage detection with estimated reference," *Mechanical Systems and Signal Processing*, vol. 164, 2022.

[J27] L. Long, M. Döhler, and S. Thöns, "Determination of structural and damage detection system influencing parameters on the value of information," *Structural Health Monitoring*, vol. 21, no. 1, 19-36, 2022.

[J28] A. Mendler, M. Döhler, and C. E. Ventura, "Sensor placement with optimal damage detectability for statistical damage detection," *Mechanical Systems and Signal Processing*, vol. 170, 2022.

[J29] A. Mendler, M. Döhler, C. E. Ventura, and L. Mevel, "Localizability of damage with statistical tests and sensitivity-based parameter clusters," *Mechanical Systems and Signal Processing*, vol. 204, 110783, 2023.

[J30] B. Merainani, B. Xiong, V. Baltazart, M. Döhler, J. Dumoulin, and Q. Zhang, "Subspace-based modal identification and uncertainty quantification from video image flows," *Journal of Sound and Vibration*, vol. 569, 117957, 2024.

[J31] A. Mendler, M. Döhler, and C. U. Grosse, "Predictive probability of detection curves based on data from undamaged structures," *Structural Health Monitoring*, accepted for publication.

Conference publications (automatic control)

[CA1] M. Döhler, X.-B. Lam, and L. Mevel, "Uncertainty quantification for stochastic subspace identification on multi-setup measurements," in *Proc. 50th IEEE Conference on Decision and Control (CDC)*, Orlando, FL, USA, 2011.

[CA2] M. Döhler and L. Mevel, "Robust subspace based fault detection," in *Proc. 18th IFAC World Congress*, Milan, Italy, 2011.

[CA3] M. Döhler and L. Mevel, "Fast multi-order stochastic subspace identification," in *Proc. 18th IFAC World Congress*, Milan, Italy, 2011.

[CA4] L. Marin, M. Döhler, D. Bernal, and L. Mevel, "Uncertainty quantification for stochastic damage localization for mechanical systems," in *Proc. 8th IFAC Symposium on Fault Detection, Diagnosis and Safety of Technical Processes (SAFE-PROCESS)*, Mexico City, Mexico, 2012.
[CA5] M. Döhler, X.-B. Lam, and L. Mevel, "Multi-order covariance computation for estimates in stochastic subspace identification using QR decompositions," in *19th IFAC World Congress*, Cape Town, South Africa, 2014.

[CA6] M. Döhler, L. Mevel, and F. Hille, "Efficient computation of minmax tests for fault isolation and their application to structural damage localization," in *Proc. 19th IFAC World Congress*, Cape Town, South Africa, 2014.

[CA7] M. Döhler and L. Mevel, "Fault isolation and quantification from Gaussian residuals with application to structural damage quantification," in *Proc. 9th IFAC Symposium on Fault Detection, Diagnosis and Safety of Technical Processes (SAFEPROCESS)*, Paris, France, 2015.

[CA8] M. Döhler, Q. Zhang, and L. Mevel, "Vibration monitoring by eigenstructure change detection based on perturbation analysis," in *Proc. 17th IFAC Symposium on System Identification (SYSID)*, Beijing, China, 2015.

[CA9] M. Döhler, Q. Zhang, and L. Mevel, "Change detection and isolation in mechanical system parameters based on perturbation analysis," in *Proc. 20th IFAC World Congress*, Toulouse, France, 2017.

[CA10] M. Bhuyan, G. Gautier, M. Döhler, Y. Lecieux, L. Mevel, and F. Schoefs, "Statistical damage localization in mechanical systems based on load vectors," in *Proc. 20th IFAC World Congress*, Toulouse, France, 2017.

[CA11] E. Viefhues, M. Döhler, F. Hille, and L. Mevel, "Asymptotic analysis of subspace-based data-driven residual for fault detection with uncertain reference," in *Proc. 10th IFAC Symposium on Fault Detection, Diagnosis and Safety of Technical Processes (SAFEPROCESS)*, Warsaw, Poland, 2018.

[CA12] A. Mendler, M. Döhler, C. E. Ventura, and L. Mevel, "Clustering of redundant parameters for fault isolation with Gaussian residuals," in *Proc. 21st IFAC World Congress*, 2020.

[CA13] E. Viefhues, M. Döhler, F. Hille, and L. Mevel, "Fault detection for linear parameter varying systems under changes in the process noise covariance," in *Proc. 21st IFAC World Congress*, 2020.

[CA14] S. Greś, M. Döhler, and L. Mevel, "Variance computation for system matrices and transfer function from input/output subspace system identification," in *Proc. 21st IFAC World Congress*, 2020.

[CA15] M. Döhler, Q. Zhang, and L. Mevel, "Damage localization in mechanical systems by Lasso regression," in *Proc. 19th IFAC Symposium on System Identification (SYSID)*, 2021.

[CA16] S. Greś, M. Döhler, and L. Mevel, "Hankel matrix-based Mahalanobis distance for fault detection robust towards changes in process noise covariance," in *Proc. 19th IFAC Symposium on System Identification (SYSID)*, 2021.
[CA17] M. Döhler, Q. Zhang, and L. Mevel, "Kalman predictor subspace residual for mechanical system damage detection," in *Proc. 11th IFAC Symposium on Fault Detection, Diagnosis and Safety of Technical Processes (SAFEPROCESS)*, Pafos, Cyprus, 2022.

[CA18] Z. Luo, B. Merainani, M. Döhler, V. Baltazart, and Q. Zhang, "High dimensional data reduction in modal analysis with stochastic subspace identification," in *Proc. 22nd IFAC World Congress*, Yokohama, Japan, 2023.

Conference publications (structural engineering/SHM)

[CS1] M. Luczak, B. Peeters, W. Szkudlarek, L. Mevel, M. Döhler, W. Ostachowicz, K. Martyniuk, and K. Branner, "Comparison of the three different approaches for damage detection in the part of the composite wind turbine blade," in *Proc. 7th International Workshop on Structural Health Monitoring (IWSHM)*, Stanford, CA, USA, 2009.

[CS2] M. Luczak, B. Peeters, M. Döhler, L. Mevel, W. Ostachowicz, P. Malinowski, T. Wandowski, and K. Branner, "Damage detection in wind turbine blade panels using three different SHM techniques," in *Proc. 28th International Modal Analysis Conference (IMAC)*, Jacksonville, FL, USA, 2010.

[CS3] M. Goursat, M. Döhler, L. Mevel, and P. Andersen, "Crystal clear SSI for operational modal analysis of aerospace vehicles," in *Proc. 28th International Modal Analysis Conference (IMAC)*, Jacksonville, FL, USA, 2010.

[CS4] D. Siegert, M. Döhler, O. Ben Mekki, L. Mevel, M. Goursat, and F. Toutlemonde, "Vibration monitoring of a small span composite bridge," in *Proc. 28th International Modal Analysis Conference (IMAC)*, Jacksonville, FL, USA, 2010.

[CS5] M. Döhler, P. Andersen, and L. Mevel, "Data merging for multi-setup operational modal analysis with data-driven SSI," in *Proc. 28th International Modal Analysis Conference (IMAC)*, Jacksonville, FL, USA, 2010.

[CS6] M. Döhler, E. Reynders, F. Magalhães, L. Mevel, G. De Roeck, and Á. Cunha, "Pre-and post-identification merging for multi-setup OMA with covariance-driven SSI," in *Proc. 28th International Modal Analysis Conference (IMAC)*, Jacksonville, FL, USA, 2010.

[CS7] M. Döhler, X.-B. Lam, and L. Mevel, "Crystal clear data fusion in subspace system identification and damage detection," in *Proc. 5th International Conference on Bridge Maintenance, Safety and Management (IABMAS)*, Philadelphia, PA, USA, 2010.

[CS8] M. Döhler, X.-B. Lam, and L. Mevel, "Confidence intervals on modal parameters in stochastic subspace identification," in *Proc. 34th International Symposium on Bridge and Structural Engineering (IABSE)*, Venice, Italy, 2010.
[CS9] M. Döhler and L. Mevel, "Modular subspace-based system identification and damage detection on large structures," in *Proc. 34th International Symposium on Bridge and Structural Engineering (IABSE)*, Venice, Italy, 2010.

[CS10] W. Zhou, H. Li, L. Mevel, M. Döhler, X.-B. Lam, C. Mao, and J. Ou, "Seismic damage assessment for a residential masonry building using aftershock monitoring of Wenchuan earthquake," in *Proc. 24th International Conference on Noise and Vibration Engineering (ISMA)*, Leuven, Belgium, 2010.

[CS11] M. Döhler, D. Siringoringo, and L. Mevel, "System identification under different ambient excitations," in *Proc. 3rd Asia-Pacific Workshop on Structural Health Monitoring (APWSHM)*, Tokyo, Japan, 2010.

[CS12] M. Döhler, B. Jaishi, L. Mevel, and J. Brownjohn, "Data fusion for system identification of the Humber Bridge," in *Proc. 29th International Modal Analysis Conference (IMAC)*, Jacksonville, FL, USA, 2011.

[CS13] M. Döhler, F. Hille, X.-B. Lam, L. Mevel, and W. Rücker, "Confidence intervals of modal parameters during progressive damage test," in *Proc. 29th International Modal Analysis Conference (IMAC)*, Jacksonville, FL, USA, 2011.

[CS14] M. Döhler, X.-B. Lam, P. Andersen, and L. Mevel, "Uncertainty quantification for stochastic subspace identification of multi-setup measurements," in *Proc. 4th International Operational Modal Analysis Conference (IOMAC)*, Istanbul, Turkey, 2011.

[CS15] M. Döhler, D. Bernal, L. Mevel, and D. Siegert, "Statistical-based model-free damage localization," in *Proc. 4th International Operational Modal Analysis Conference (IOMAC)*, Istanbul, Turkey, 2011.

[CS16] M. Döhler, L. Mevel, and D. Bernal, "Statistical subspace-based model-free damage localization," in *Proc. Engineering Mechanics Institute Conference (EMI)*, Boston, MA, USA, 2011.

[CS17] F. Hille, M. Döhler, L. Mevel, and W. Rücker, "Subspace based damage detection methods on a prestressed concrete bridge," in *Proc. 8th International Conference on Structural Dynamics (EURODYN)*, Leuven, Belgium, 2011.

[CS18] M. Döhler, F. Hille, L. Mevel, and W. Rücker, "Structural health monitoring during progressive damage test of S101 Bridge," in *Proc. 8th International Workshop on Structural Health Monitoring (IWSHM)*, Stanford, CA, USA, 2011.

[CS19] J. Strout, M. Döhler, D. Bernal, and L. Mevel, "Changes in the statistics of ambient excitations in the performance of two damage detection schemes," in *Proc. 30th International Modal Analysis Conference (IMAC)*, Jacksonville, FL, USA, 2012.

[CS20] M. Döhler, P. Andersen, and L. Mevel, "Operational modal analysis using a fast stochastic subspace identification method," in *Proc. 30th International Modal Analysis Conference (IMAC)*, Jacksonville, FL, USA, 2012.

[CS21] M. Döhler, L. Mevel, and D. Siegert, "Statistical subspace-based damage detection under changing ambient excitation," in *Proc. 5th European Conference on Structural Control (EACS)*, Genoa, Italy, 2012.
[CS22] P. Mellinger, M. Döhler, L. Mevel, and G. Broux, "Data fusion for frequency domain stochastic subspace identification," in *Proc. 5th European Conference on Structural Control (EACS)*, Genoa, Italy, 2012.

[CS23] S. Thöns and M. Döhler, "Structural reliability updating with stochastic subspace damage detection information," in *Proc. 5th European Conference on Structural Control (EACS)*, Genoa, Italy, 2012.

[CS24] M. Döhler, L. Mevel, and P. Andersen, "Efficient uncertainty computation for modal parameters in stochastic subspace identification," in *Proc. 4th International Conference on Uncertainty in Structural Dynamics (USD)*, Leuven, Belgium, 2012.

[CS25] D. Bernal, S. Mozaffari Kojidi, K. Kwan, and M. Döhler, "Damping identification in buildings from earthquake records," in *Proc. CSMIP Annual Seminar on Utilization of Strong Motion Data*, Sacramento, CA, USA, 2012.

[CS26] L. Marin, M. Döhler, D. Bernal, and L. Mevel, "Damage localization using a statistical test on residuals from the SDDLV approach," in *Proc. 31st International Modal Analysis Conference (IMAC)*, Garden Grove, CA, USA, 2013.

[CS27] S. Mozaffari Kojidi, M. Döhler, D. Bernal, and Y. Liu, "Linear projection techniques in damage detection under a changing environment," in *Proc. 31st International Modal Analysis Conference (IMAC)*, Garden Grove, CA, USA, 2013.

[CS28] M. Döhler, K. Kwan, and D. Bernal, "Optimal sensor placement with a statistical criterion for subspace-based damage detection," in *Proc. 31st International Modal Analysis Conference (IMAC)*, Garden Grove, CA, USA, 2013.

[CS29] S. Mozaffari Kojidi, D. Bernal, and M. Döhler, "Change detection under environmental variability," in *Proc. 5th International Operational Modal Analysis Conference (IOMAC)*, Guimarães, Portugal, 2013.

[CS30] M. Döhler, L. Marin, D. Bernal, and L. Mevel, "Comparison of two statistical damage localization approaches," in *Proc. 5th International Operational Modal Analysis Conference (IOMAC)*, Guimarães, Portugal, 2013.

[CS31] L. Marin, M. Döhler, D. Bernal, and L. Mevel, "Statistical based decision making for damage localization with influence lines," in *Proc. 9th International Workshop on Structural Health Monitoring (IWSHM)*, Stanford, CA, USA, 2013.

[CS32] V. Le Cam, M. Döhler, M. Le Pen, and L. Mevel, "Embedded modal analysis algorithms on the smart wireless sensor platform PEGASE," in *Proc. 9th International Workshop on Structural Health Monitoring (IWSHM)*, Stanford, CA, USA, 2013.

[CS33] M. Döhler and F. Hille, "Subspace-based damage detection on steel frame structure under changing excitation," in *Proc. 32nd International Modal Analysis Conference (IMAC)*, Orlando, FL, USA, 2014.

[CS34] D. Bernal, M. Döhler, and D. Parker, "Examination of two sensor placements schemes in damage detection," in *Proc. 7th European Workshop on Structural Health Monitoring (EWSHM)*, Nantes, France, 2014.
[CS35] C. Ventura, P. Andersen, L. Mevel, and M. Döhler, "Structural Health Monitoring of the Pitt River Bridge in British Columbia, Canada," in *Proc. 6th World Conference on Structural Control and Monitoring*, Barcelona, Spain, 2014.

[CS36] M. Döhler, L. Marin, L. Mevel, and D. Bernal, "Operational modal analysis with uncertainty quantification for SDDLV-based damage localization," in *Proc. 4th French congress Analyse Vibratoire Expérimentale*, Blois, France, 2014.

[CS37] S. Allahdadian, C. Ventura, P. Andersen, L. Mevel, and M. Döhler, "Sensitivity evaluation of subspace-based damage detection method to different types of damage," in *Proc. 33rd International Modal Analysis Conference (IMAC)*, Orlando, FL, USA, 2015.

[CS38] S. Allahdadian, C. Ventura, P. Andersen, L. Mevel, and M. Döhler, "Investigation on the sensitivity of subspace based damage detection technique to damage and noise levels," in *Proc. 6th International Operational Modal Analysis Conference (IOMAC)*, Gijón, Spain, 2015.

[CS39] S. Allahdadian, C. Ventura, P. Andersen, L. Mevel, and M. Döhler, "Subspace based damage detection technique: Investigation on the effect of number of samples," in *Proc. 11th Canadian Conference on Earthquake Engineering*, Victoria, BC, Canada, 2015.

[CS40] S. Allahdadian, M. Döhler, C. Ventura, and L. Mevel, "On the influence of sample length and measurement noise on the stochastic subspace damage detection technique," in *Proc. 34th International Modal Analysis Conference (IMAC)*, Orlando, FL, USA, 2016.

[CS41] M. Döhler and S. Thöns, "Efficient structural system reliability updating with subspace-based damage detection information," in *Proc. 8th European Workshop on Structural Health Monitoring (EWSHM)*, Bilbao, Spain, 2016.

[CS42] G. Gautier, L. Mevel, J.-M. Mencik, M. Döhler, and R. Serra, "Parameter uncertainties quantification for finite element based subspace fitting approaches," in *Proc. 8th European Workshop on Structural Health Monitoring (EWSHM)*, Bilbao, Spain, 2016.

[CS43] V. Le Cam, M. Döhler, M. Le Pen, I. Guéguen, and L. Mevel, "Embedded subspace-based modal analysis and uncertainty quantification on wireless sensor platform PEGASE," in *Proc. 8th European Workshop on Structural Health Monitoring (EWSHM)*, Bilbao, Spain, 2016.

[CS44] M. Bhuyan, M. Döhler, and L. Mevel, "Statistical damage localization with stochastic load vectors using multiple mode sets," in *Proc. 8th European Workshop on Structural Health Monitoring (EWSHM)*, Bilbao, Spain, 2016.

[CS45] M. Limongelli, E. Chatzi, M. Döhler, G. Lombaert, and E. Reynders, "Towards extraction of vibration-based damage indicators," in *Proc. 8th European Workshop on Structural Health Monitoring (EWSHM)*, Bilbao, Spain, 2016.
[CS46] P. Mellinger, M. Döhler, and L. Mevel, "Variance estimation of modal parameters from input/output covariance-driven subspace identification," in *Proc. 27th Conference on Noise and Vibration Engineering (ISMA)*, Leuven, Belgium, 2016.

[CS47] M. Bhuyan, E. Viefhues, M. Döhler, Y. Lecieux, L. Mevel, F. Hille, and F. Schoefs, "Output-only subspace and transfer matrix-based damage localization and quantification," in *Proc. 7th International Operational Modal Analysis Conference (IOMAC)*, Ingolstadt, Germany, 2017.

[CS48] M. Döhler, P. Andersen, and L. Mevel, "Variance computation of modal parameter estimates from UPC subspace identification," in *Proc. 7th International Operational Modal Analysis Conference (IOMAC)*, Ingolstadt, Germany, 2017.

[CS49] D. Straub, E. Chatzi, E. Bismut, W. Courage, M. Döhler, M. Faber, J. Köhler, G. Lombaert, P. Omenzetter, M. Pozzi, S. Thöns, D. Val, H. Wenzel, and D. Zonta, "Value of information: A roadmap to quantifying the benefit of structural health monitoring," in *Proc. 12th International Conference on Structural Safety & Reliability (ICOSSAR)*, 2017.

[CS50] S. Greś, M. D. Ulriksen, M. Döhler, R. J. Johansen, P. Andersen, L. Damkilde, and S. A. Nielsen, "Statistical methods for damage detection applied to civil structures," in *Proc. 10th International Conference on Structural Dynamics (EURODYN)*, Rome, Italy, 2017.

[CS51] E. Viefhues, M. Döhler, F. Hille, and L. Mevel, "Stochastic subspace-based damage detection with uncertainty in the reference null space," in *Proc. 11th International Workshop on Structural Health Monitoring (IWSHM)*, Stanford, CA, USA, 2017.

[CS52] M. Bhuyan, S. Allahdadian, M. Döhler, Y. Lecieux, L. Mevel, F. Schoefs, and C. Ventura, "Transfer matrices-based output-only statistical damage localization and quantification on a frame structure," in *Proc. 11th International Workshop on Structural Health Monitoring (IWSHM)*, Stanford, CA, USA, 2017.

[CS53] L. Long, S. Thöns, and M. Döhler, "Damage detection and deteriorating structural systems," in *Proc. 11th International Workshop on Structural Health Monitoring (IWSHM)*, Stanford, CA, USA, 2017.

[CS54] S. Thöns, M. P. Limongelli, A. M. Ivankovic, M. Faber, D. Val, M. Chryssanthopoulos, G. Lombaert, M. Döhler, D. Straub, E. Chatzi, J. Köhler, H. Wenzel, and J. Sorensen, "Progress of the COST Action TU1402 on the quantification of the value of structural health monitoring," in *Proc. 11th International Workshop on Structural Health Monitoring (IWSHM)*, Stanford, CA, USA, 2017.

[CS55] S. Allahdadian, M. Döhler, C. Ventura, and L. Mevel, "Damage localization of a real structure using the statistical subspace damage localization method," in *Proc. 11th International Workshop on Structural Health Monitoring*, Stanford, CA, USA, 2017.

[CS56] L. Long, S. Thöns, and M. Döhler, "The effects of SHM system parameters on the value of damage detection information," in *Proc. 9th European Workshop on Structural Health Monitoring (EWSHM)*, Manchester, UK, 2018.