Abstract

This study explores a change detection method in modal properties to automate and generalize in-service damage detection for vibration-based structural health monitoring of bridges. The noisy conditions caused by ambient loading pose difficulty for in-service damage detection because the load-induced noise often masks the difference in the modal properties. The proposed method directly converts measured time series into a simplified anomaly indicator robust against load-induced noise. This study adopts a vector autoregressive model to represent the vibration of bridges. Bayesian inference produces a posterior probability distribution function of the model parameters. Principal component analysis extracts a subspace comparable to the modal properties in the model parameters. Bayesian hypothesis testing quantifies anomalies in the extracted subspace. The feasibility of the proposed method is assessed with vibration data from field experiments conducted on an actual steel truss bridge. The field experiment includes damage severing the truss members. The modal frequencies and mode shapes estimated from the principal component analysis correspond well to earlier reported results. The proposed damage detection method successfully indicated all damage considered in the experiment.

1. Introduction

1.1. Background

Aging infrastructure management is a crucially important issue confronted by civil engineering professionals. Today, they adopt periodic visual inspection-based bridge management as a fundamental method of investigating the structural integrity of bridges. Structural health monitoring (SHM) provides a research field for nondestructive evaluation based on physical measurements and computer analyses to complement the existing inspection methods. Especially, vibration-based SHM has the advantage that a few sensors enable assessment of the condition of the whole structure because a local stiffness change affects the global dynamic characteristics represented by modal properties [1–3]. Two commonly used means to identify modal properties are input-output estimation and output-only estimation [3]. For input-output estimation, bridges are excited under controlled input signals. Contrary, bridges are excited by uncontrolled ambient loads such as traffic, wind, and ground motion for output-only estimation. We also refer to the output-only methods as operational modal analysis (OMA) methods because they enable the identification of modal properties solely from the measurable outputs under operation without traffic interruption. The ambient loads are usually modeled as steady white noise in most OMA methods since the external force is unknown in advance. We expect to apply OMA for damage detection as a screening method for numerous bridges without traffic interruption. However, as discussed below, noisy conditions caused by ambient loads still require effort for calibration and computation to automate damage detection.

Aiming damage detection using OMA, most earlier studies have adopted statistic change detection in modal properties, such as modal frequency, mode shapes, and damping ratios [4–14]. Among them, Gudmundson [6] theoretically demonstrated that cracks and notches in a beam alter modal frequencies. Pandey et al. [7, 8] proposed to detect local changes in flexural stiffness using the curvature of mode shapes. Deramaeker et al. [9] demonstrated the feasibility of output-only modal estimation for vibration-based SHM under environmental changes with sensitivity analysis using a numerical model. For actual bridges in the field, the authors [10, 11] reported that they detected changes in modal properties caused by damage from vehicle-induced vibration. Several studies investigated the feasibility of OMA for the long-term monitoring of bridges under operation. Peeters and De Roeck [12] applied black-box models for one-year bridge monitoring data and distinguished damage-induced frequency changes from environmental-induced changes. Ni et al. [13] reported a change in dynamic behavior between ambient vibration and typhoon-induced vibration. Although environmental fluctuation changes modal properties, damage-induced fluctuation would be distinguished from environmental effects if we have enough prior knowledge and apply an appropriate change detection procedure. However, engineers must also consider that the load-induced noise often masks modal estimation because it engenders unrealistic modal properties estimators [10, 11, 14]. That is because ambient loads are modeled as steady white noise in most OMA methods, whereas ambient loads are often unsteady and colored noise. Therefore, they need to remove the unrealistic estimators caused by the masking effect of ambient loads.

The following feature extraction techniques are required to avoid the masking effect and distinguish the physically meaningful modes from the meaningless estimators: the frequency domain OMA requires preliminary calibration for peak-pick [15]. In time domain OMA, they usually need thresholds of modal validation criteria for stabilization diagrams [3]. In recent studies, Mao et al. [16] proposed to automate modal identification with hierarchical clustering and principal component analysis (PCA) in addition to the classical stabilization diagrams. In the previous study, they applied PCA to reduce the dimension of modal validation criteria beforehand the clustering process. Engineers can compare modal properties estimated with the abovementioned feature extraction techniques for damage detection. However, they still need to apply statistical methodology such as hypothesis testing because each modal property involves variance caused by the estimation error. Also, they are often required to choose modal properties for comparison because the modal properties irrelevant to the damage decrease the sensitivity of the hypothesis test. The comparison requires engineers’ experience both for hypothesis testing and for modal analysis, and thus the results of damage detection depend on the skill of the engineers. A subjective judgment made by inexperienced engineers may cause a false alert or overlook failure. Accordingly, it is still challenging to generalize and automate the procedure required for change detection in modal properties.

Nair et al. [17] proposed omitting modal estimation procedures for damage detection to automate damage detection. They applied anomaly indicators converted from a univariate AR model instead of modal properties. The feasibility of the anomaly indicator has been examined for a modeled bridge [18] and an actual steel truss bridge [19]. Zhang [20] also proposed a statistical diagnosis procedure by change detection in the AR model. The proposed procedure examines the difference in residuals in the AR models. For multipoint measurement, the authors [21] expanded these previous studies and proposed an anomaly indicator calculated from a vector autoregressive (VAR) model instead of a univariate AR model. In the study, the authors carefully chose the AR order to avoid the masking effect caused by the load-induced noise. That is because the load-induced noise produces physically meaningless estimators in the AR coefficients. The practical application of AR-model-based anomaly indicators is consequently limited so far. For an anomaly indicator robust against ambient loads, we desire to extract physically meaningful parameters from the VAR model and cancel out the masking effect caused by the meaningless VAR coefficients. This study applies Bayesian statistics to quantify the uncertainty involved in VAR coefficients and proposes a methodology to cancel them out. In a Bayesian manner, the procedure to adjust the AR order is also automated and generalized. Furthermore, this paper applies PCA to extract the physically meaningful parameters relevant to classical modal analysis theory.

Bayesian statistics is a widely used methodology to quantify the uncertainty in stochastic model parameters. Several previous studies applied Bayesian inference [22] to quantify the uncertainty of modal properties. For instance, Yuen and Katafygiotis [23] and Au et al. [24] proposed a Bayesian modal analysis method. In the previous study, the posterior distribution of the parameters represents their uncertainty according to observed data. Lam et al. [25] applied this method to model updating, and the updated model was subsequently adopted for model-based damage detection. Some other previous studies applied Bayesian methodologies in various approaches for damage detection. For instance, Jiang and Mahadevan [26] proposed a Bayesian hypothesis testing [27] adopting a dynamic fuzzy wavelet neural network method for nonparametric damage detection. Sun et al. [28] proposed damage detection using a Bayesian information model, including local information such as static strain and deflection. Figueiredo et al. [29] applied Mahalanobis distance obtained from Markov chain Monte Carlo (MCMC) techniques as an anomaly indicator. Dzunic et al. [30] applied a Bayesian state-space approach and investigated multiclass classification for damage detection. The Bayesian stochastics is also available for machine learning providing decision-making strategies for damage detection. For instance, Silva et al. [31] adopted neural networks with autoencoders to extract damage-sensitive features. Arangio and Beck [32] utilized Bayesian neural networks and discussed their feasibility for SHM. In the recent advances in uncertainty quantification and machine learning, stochastic models tend to be detailed and complicated. Accordingly, some of the parameters that consist of the models and threshold are often required to be predefined. Also, most cases for recent Bayesian applications adopt iterative optimization methods such as MCMC techniques. The computational resources needed in optimization are considerable for prompt decision-making in SHM. As discussed later, Bayesian inference for the VAR model applied in this paper enables much simpler computation because it has posterior distribution in closed form solution.

Recently, the novelty detection approach in machine learning has been applied to structural damage detection [33, 34]. In this sense, one-class classification formulates structural damage detection without prior knowledge of damages. This study concerns parametric change detection utilizing this approach. Damage detection methods proposed by Figueiredo et al. [29] and Dzunic et al. [30] mentioned above also apply this approach. Dimension reduction by PCA is one of the popular methodologies to extract physically meaningful parameters for change detection [35–37]. For instance, Akintude et al. [35] demonstrated that proper orthogonal modes from PCA of strain time history indicated damage in a full-scale bridge deck mock-up using variable vehicle loads and speeds. In the previous study by Mao et al. [16] mentioned above, PCA is applied to extract physically meaningful modes from stabilization diagrams. For anomaly detection of modal properties, Ozdagli and Kousoukos [36] investigated an indicator to detect changes in modal properties using reduced dimensions by PCA. PCA was also applied to distinguish damage-induced frequency change from environment-induced change. For instance, Sen et al. [37] applied PCA to modal frequencies to decouple structural damage and environmental effects. Furthermore, a classical modal analysis theory shows that the modal properties are also related to the PCA for a vector of the measured time series, as discussed in Appendix A. This study applies this result to extract damage-sensitive parameters for anomaly detection.

1.2. Outline of the Proposed Method

This study aims to provide an easy-to-use framework for change detection in modal properties with minimum adjustment. This study proposes a change detection method so that it meets all of the following conditions:(i)An anomaly indicator is formulated based on hypothesis testing for convincing change detection(ii)To automate and generalize the procedure for damage detection, the anomaly indicator is directly converted from the measured time series(iii)Bayesian statistics is applied to quantify the uncertainty involved in the stochastic model(iv)For fast computation, no iterative method is adopted for Bayesian inference(v)Damage-sensitive parameters comparable to modal properties are extracted to avoid errors caused by parameters unrelated to modal properties

This study adopts a VAR model as a bridge vibration model under ambient loads. Even though a VAR model merely provides a rough approximation of the structural responses, it still has an advantage in simple and fast computation. As discussed in Section 1.1, unsteady load-induced noise engenders meaningless estimators in the VAR coefficient. This study proposes the following framework to reduce the masking effect caused by ambient loads.

Bayesian inference for the VAR coefficients and the covariance matrix of the error term provide posterior distribution of the parameters. Because the posterior distribution is in closed form solution, the Bayesian inference is available without the iterative optimization. In the Bayesian inference, the prior distribution represents noninformative initial knowledge before measurement and the posterior distribution represents updated knowledge after measurement. Hyperparameters calculated from the measured time series represent the posterior distribution of the VAR coefficients and the covariance matrix. The posterior distribution quantifies the most probable values and uncertainty for each regressive parameter from the observation. An orthogonal basis derived from the PCA transforms the posterior distribution of VAR coefficients. The damage-sensitive parameters appear in a parametric subspace whose dimension the PCA reduces. Appendix A shows that the extracted damage-sensitive parameters are comparable to the modal properties estimated in a classical modal analysis method. Accordingly, the transformed posterior distribution separates the VAR coefficients into damage-sensitive and meaningless parameters. Bayesian hypothesis testing [27] is applied to detect changes in the damage-sensitive parameters. The hypothesis testing adopts two hypotheses representing the healthy and damaged state of the bridge. For both hypotheses, evidence functions represent the likelihood of the observed data. Here, one evidence function presumes changes in the damage-sensitive parameters, and the other assumes no change. Each evidence function marginalizes the uncertainty of the other meaningless parameters with integration. A Bayes factor, a ratio of the two evidence functions, indicates the changes. Accordingly, the proposed Bayes factor reduces the masking effect caused by ambient loads. In this study, the threshold for the log-scaled Bayes factor to detect damage is fixed at 0, as discussed in Section 3.

The above procedure comprises the following two steps: Step 1: Reference Modeling Acquire acceleration time series from a bridge under healthy conditions (designated as the reference dataset). According to the Bayesian inference, transform the reference dataset into hyperparameters. Find the optimal AR order by the Bayesian information criterion (BIC) [22, 38]. Extract the damage-sensitive parameters by PCA. Step 2: Damage Detection Acquire a new acceleration time series from the same bridge under unknown damage conditions (designated as the test dataset). Calculate hyperparameters of the evidence functions for the two hypotheses representing healthy and damaged states. The Bayes factor exceeds the threshold if it detects any changes in the damage-sensitive parameters between the reference and test datasets.

Engineers can convert damage-sensitive parameters to modal properties through a relevant state-space model [39, 40] just after Step 1. This procedure helps them interpret the physical meaning of the damage-sensitive parameters. They can omit this procedure from the framework if they do not need the physical interpretation.

Figure 1 presents an outline of the proposed damage detection method. Here, the equation numbers appear in Section 2 and Section 3. Section 2 describes optimizing the stochastic model for Step 1 and modal estimation. In Section 3, we formulate hypothesis testing for Step 2. We validate the proposed method in Section 4 with a case study of a steel truss bridge in the field [41]. Furthermore, Appendix B provides a case study using an existing approach, i.e., Bayesian hypothesis testing without principal component analysis, for comparison.

Figure 1 

Outline of the proposed method.

2. Reference Modeling

2.1. VAR Model Expansion

This section presents a general description of the VAR model expanded to Bayesian statistics. Let a synchronized time series of the acceleration be measured at m measurement locations on a bridge. A zero-mean time series is produced by subtracting the mean value from the measured acceleration. Let () denote a column vector at the k-th time step of the zero-mean time series with components corresponding to the measurement locations. As mentioned in Section 1, most of the OMA methods assume the ambient loads as steady white noise. The VAR model with sufficient numbers of model order approximates the time series obtained from a linear structural dynamic system excited by white noise [42]. The VAR model is given as

In that equation, represents the i-th AR coefficient matrix, and denotes a zero-mean Gaussian white noise vector. Letting and , equation (1) leads to the following probability density function (PDF).

Therein, denotes the covariance matrix of ; denotes the PDF of ; denotes the conditional PDF of conditioned by and ; denotes the PDF of following the multivariate Gaussian distribution with expectation and covariance matrix . The conditional PDF in equation (2) is the likelihood function of for parameters and . According to the Bayes theorem, the PDF of and conditioned by is given aswhere and represent the prior and posterior PDF of and , respectively; is a constant given as the marginal likelihood of .

For Bayesian inference, this study adopts conjugate priors [22], which produce posterior distributions with invariant functional forms. The conjugate priors have an advantage for fast computation because they do not require MCMC. For a linear regressive model, the conjugate prior for and is given as the following matrix-normal inverse-Wishart distribution [30].

Therein, and stand, respectively, for the matrix-normal and inverse-Wishart distributions aswhere represents the trace of a matrix, denotes the determinant of a matrix, and is the multivariate gamma function. , , , and are hyperparameters, i.e., the parameters which determine the functional form of the PDF for the parameters and . Equations (2)–(4) engender the following posterior distribution [30].

In equation (6), each hyperparameter is given as follows:

Equations (2)–(4) and (6) lead the following marginal likelihood.

Herein, represents that the right-hand side is proportional to the left-hand side.

2.2. Bayesian Inference for the Reference Dataset

Let sets of time series of acceleration be acquired as the reference dataset: stands for each of the time series where is its data length. Let represent the reference dataset as , where the subscript “r” denotes “reference.” Before observing the reference dataset, we do not have any information about the value of the model parameters. Therefore, it is reasonable to adopt a so-called noninformative prior [22, 43], which represents unknown initial knowledge about the parameters. In equations (7)–(10), this study applies hyperparameters to represent the noninformative prior. As far as the abovementioned hyperparameters are adopted, the value of does not alter the following calculation. The proposed noninformative prior and equations (7)–(10) generate the following posterior distribution:

In equation (12), each hyperparameter is calculated aswhere . The posterior distribution in equation (12) represents the updated knowledge after measuring the reference dataset and quantifies the uncertainty for each regressive parameter based on the observation. The time series in the reference dataset directly produce the hyperparameters in equations (13)–(16). The summations in equations (13), (14), and (16) are approximated with the following autocorrelation for fast computation:

Herein, represents the total data length of the reference dataset, i.e., . Equations (13), (14), and (16) are approximated using autocorrelation functions, as shown in the following equations:

To find appropriate AR order p, the following BIC [22, 38] is adopted.

Therein, and are the most likelihood estimators (MLEs) of and for the reference dataset. The most likelihood estimators can be rewritten using the hyperparameters as

According to the well-known nature of the Gaussian model, in equation (19) is given as presented in the following equation:

Using the hyperparameters and ignoring the constant, BIC in equation (19) is rewritten as

The AR order providing the lower BIC represents the better regressive model representing the reference dataset. Accordingly, the optimal AR order is automatically calculated from the reference dataset.

2.3. Principal Component Analysis

The posterior distribution for is transformed into independent parameters for further analysis. For the PDF given in equation (12), the orthogonal basis for the transformation is obtained from the eigenvalue decomposition of the hyperparameter aswhere is the diagonal matrix consisting of the eigenvalues; is the orthogonal matrix consisting of the eigenvectors; stands for the diagonal matrix consisting of the diagonal elements . and stand, respectively, for the i-th column of and the i-th diagonal member of . For simplicity, let the eigenvalues be in descending order, i.e., . The orthogonal basis transforms the posterior distribution for as follows: Letting and , the PDF in equation (12) is transformed aswhere and . Equation (25) shows that each vector is independent of the others.

The transformation described previously is comparable to the PCA of the vector . That is because equation (13) relates the abovementioned orthogonal basis and the covariance matrix of the vector as follows:

This equation demonstrates that the vectors are the eigenvectors of the covariance matrix of . Therefore, the i-th principal components of are given as . According to the PCA theory, vector is approximated with the first several principal components. The following contribution ratio is used widely to evaluate how principal components contribute to the approximation.

The abovementioned contribution ratio is applied to distinguish the damage-sensitive parameters from the other parameters.

The following steady state-space model is examined to clarify the relationship between modal properties and the PCA.

Here, is a state vector representing structural response, and matrices , , and are constant matrices. The modal frequencies, damping ratios, and mode shapes are estimated from matrices and . As discussed in Appendix A, the eigensystem realization algorithm (ERA) [40] produces the estimators of and as

In those equations, represents a vector that consists of the a-th to b-th elements of .

For the reference dataset, the modal properties are estimated by substituting the MLE given in equation (20) to . The MLE of matrix is given as

Modal properties are estimated with eigenvalue analysis of matrix aswhere is one of the eigenvalues of and is the corresponding eigenvector. The angular frequency , damping ratio , and the mode shape s are given as follows:where , , and represent the absolute value, real part, and phase angle, respectively; stands for the sampling interval.

Equations (30) and (31) show that the AR coefficients relate with the modal properties only on the vector space consisting of the orthogonal basis . Namely, only q vectors … represent the modal properties, and the other vectors are less relevant to them. The orthogonal transformation consequently separates physically meaningful and meaningless parameters. This study refers to the vectors … as the damage-sensitive parameters. Let the damage-sensitive parameters and the other parameters be gathered in matrices, respectively, as and . Equation (25) produces the following PDFs for , , and .

3. Damage Detection

3.1. Hypothesis Testing

Bayesian hypothesis testing provides a framework to compare two competing statistical models representing the observed data [27]. This study adopts two statistical models, respectively, representing the healthy and damaged bridge. Following the hypothesis testing precedent, the healthy and damaged bridge models are designated as the “null hypothesis” and the “alternative hypothesis.”

Let denote the test dataset for the hypothesis testing. Herein, subscript “t” stands for “test.” Let and represent the damage-sensitive parameters for the reference and test datasets, respectively. This study formulates the hypotheses according to the following conditions:(i)The damage-sensitive parameters for the null hypothesis are consistent between the reference and test datasets, i.e., (ii)The damage-sensitive parameters for the alternative hypothesis differ between the reference and test datasets, i.e., (iii)Parameters and for both hypotheses are consistent between the reference and test datasets

Under the null hypothesis, no parameter is altered between the reference and the test datasets. Therefore, the PDF of the parameters for the null hypothesis is identical to equation (12), i.e.,

Under the alternative hypothesis, damage-sensitive parameters are altered from the reference dataset and unknown before observing the test dataset. Bayesian statistics represents the parameters with less knowledge as PDFs with broad variances. As shown in equation (37), each vector has variance in the posterior distribution. The eigenvalues of are usually large enough, i.e., (i = 1, …, q) because sufficient data length is available for the reference dataset. So, the alternative hypothesis adopts the following PDF with broad variances instead of equation (37):

Here, 0 is adopted as the mean value of the damage-sensitive parameters to avoid bias derived from the reference dataset.

Consequently, the PDF for parameters and are given aswhere

Letting , , and for simplicity, the PDFs for the null and alternative hypotheses are presented as follows:

Using the PDF in equation (44), evidence functions of for the hypothesis testing are

The only difference between the two hypotheses is the distribution of the damage-sensitive parameters relevant to the modal properties. The uncertainty involved in the other parameters is marginalized and canceled out by the integration in equation (45).

The Bayes factor is the ratio of the two evidence functions shown as follows:

According to the well-known scale as the likelihood ratio test statistics, Kass and Raftery [27] interpreted the Bayes factor on twice the natural logarithm scale as

When the test dataset is in favor of the alternative hypothesis, its evidence function is higher against the null hypothesis, i.e., . So, if the damage-sensitive parameters differ from the reference dataset, likely exceeds 0, and vice versa. Therefore, this study adopts as the threshold for damage detection where the model evidence for the null and the alternative hypotheses are equivalent, i.e., . Kass and Raftery [27] proposed the scale shown in Table 1. For example, if is greater than 10, the evidence against the null hypothesis is interpreted as “very strong.”

3.2. Global Damage Indicator

Let represent a time series observed for hypothesis testing, where stands for the data length. This section adopts the whole time series as the test dataset, i.e., . Equation (11) gives the following evidence function:

In equation (48), , , and stand for the hyperparameters of the posterior distribution after test dataset observation, i.e., the hyperparameters of . The following equations calculate each hyperparameter:

Here, . Accordingly, equation (48) gives the following log-scaled Bayes factor.