Bayesian Prediction of Pre-Stressed Concrete Bridge Deflection Using Finite Element Analysis

Vertical deflection has been emphasized as an important safety indicator in the management of railway bridges. Therefore, various standards and studies have suggested physics-based models for predicting the time-dependent deflection of railway bridges. However, these approaches may be limited by model errors caused by uncertainties in various factors, such as material properties, creep coefficient, and temperature. This study proposes a new Bayesian method that employs both a finite element model and actual measurement data. To overcome the limitations of an imperfect finite element model and a shortage of data, Gaussian process regression is introduced and modified to consider both, the finite element analysis results and actual measurement data. In addition, the probabilistic prediction model can be updated whenever additional measurement data is available. In this manner, a probabilistic prediction model, that is customized to a target bridge, can be obtained. The proposed method is applied to a pre-stressed concrete railway bridge in the construction stage in the Republic of Korea, as an example of a bridge for which accurate time-dependent deflection is difficult to predict, and measurement data are insufficient. Probabilistic prediction models are successfully derived by applying the proposed method, and the corresponding prediction results agree with the actual measurements, even though the bridge experienced large downward deflections during the construction stage. In addition, the practical uses of the prediction models are discussed.


Introduction
The vertical deflection of railway bridges is an important safety indicator that is used in the safety management of structures and trains [1,2]. In particular, high-speed trains can be derailed owing to railway defects that can be caused by excessive bridge deflection. Therefore, the monitoring and management of the vertical deflection of railway bridges have been emphasized in several standards, such as the Design Guide for Steel Railway Bridges of the UK (2004), UIC CODE 518 OR (2009), and the Guideline of Track Maintenance of Korea Rail Network Authority (2016), which suggest monitoring for the acceptable vertical deflection at mid-span [3][4][5]. For example, the Korean railway bridge guideline specifies allowable maximum deflections with respect to the span length (L) and train velocity. For instance, for train speeds between 200 km/h and 270 km/h, the allowable threshold of vertical deflection with a span length is shorter than 20 m is L/1300 m [5].
Pre-stressed concrete (PSC) girders have been widely adopted as superstructures for high-speed railway bridges; the long-term deflection of bridges with PSC girders is known to be dependent on D = (X, y) = (x ij , y i ) i = 1, . . . , N D ; j = 1, . . . , N x (1) where x ij is the element of the training input matrix X, whose matrix size is N D by N x , and y i is the element in the training output vector y, whose size is N D by 1. The input-output relationship is expressed as, where φ(·) is the vector of N dim projection functions (e.g., φ(z) = [1, z, z 2 , z 3 ] for the three-dimensional polynomial regression of a scalar input z or φ(z) = [1, sin(2πz), cos(2πz)] for the first-order trigonometric regression), f (X) is the regression function of X, which can be calculated as the product of the projected input matrix φ(X), and the vector of weights w (whose sizes are N D by N dim and N dim by 1, respectively), and ε is the vector of observational errors (whose vector size is N D by 1). Here, the observational errors are assumed to follow an independent and identically distributed Gaussian distribution with a zero-mean and constant variance of the noise, σ n 2 for all inputs (i.e., a homoscedastic noise assumption) as: To infer the weight vector w based on actual measurements in the Bayesian view, the prior assumption on the weight vector is set to be a zero-mean Gaussian with the covariance matrix Σ p : which enables the prior prediction on the mean µ X* and the standard deviation σ X* based on the given input matrix X * as: From the independent and identically distributed Gaussian process assumption of observations having homoscedastic noise, the likelihood function can be constructed as the multiplication of normal probability density functions (PDFs) as follows, where y i is an i-th element of the output vector, x i is the corresponding input vector, · 2 is an L2 norm which denotes the Euclidean distance. Consequently, Equation (6) indicates the following relation: Meanwhile, from the Bayes' rule, the posterior probability on the weight vector is expressed as: p w X, y = p( X,y|w)p(w) p(X,y) = p(y|X,w )p(X,y) p(y|X) It can be seen from Equation (8) that the posterior probability p(w|X,y) is proportional to the prior probability p(w) and the likelihood p(y|X,w), which can be expressed as follows: p w X, y ∝ p(y|X, w )p(w). (9) From Equations (4) and (6), the posterior probability in Equation (9) is seen to be proportional to the multiplication of two multi-variate Gaussian PDFs: That is [31], and For the given input matrix X * , the predictive mean φ(X * )·w is [31], and the predictive variance is [31], which is based on the following matric inversion [33]: In summary, the predictive mean and variance with the posterior inference on the weights vector are: Here, the term φ(·)Σ p φ(·) T is the covariance matrix of the projected inputs, and when expressed with the capital letter "K," Equation (17) is arranged as:

Covariance Matrix Design Using Kernel Functions
In GPR, the covariance matrix of the projected inputs is numerically estimated using kernel functions. That is, the covariance of two input vectors X a and X b is defined by the multiplication of two standard deviations and one correlation coefficient, and is expressed with a kernel function k as follows, where σ x a and σ x b are the standard deviations of X a , and X b , respectively, and ρ x a x b is the correlation coefficient of X a and X b . The covariance matrix, based on kernel functions, enables the GP model to fit data even though the data exhibits nonlinearity [31]. If data have the same variance, then the input data are homoscedastic, and the standard deviation must be a constant that is the square root of the homoscedastic variance. Meanwhile, a correlation coefficient between the two stationary data is generally expressed to be inversely proportional to the Euclidean distance for stationary data, and one representative is the squared exponential kernel (k SE ) [31], where θ 1 2 is the variance hyperparameter that controls the vertical scale of the function change, and θ 2 2 is the length-scale hyperparameter that controls the horizontal scale of the function change [31].
The variance of noise (σ n 2 ) in Equation (2) is often expressed by a kernel function (k var ) using the Kronecker delta function, where θ 3 2 is the variance hyperparameter for the input noise. In addition, a bias is introduced using a constant kernel (k c ) with the hyperparameter θ 4 : In this study, a new kernel is suggested by combining the following kernels: These hyperparameters strongly affect the performance of GPR, and their optimal values are estimated based on the given observation data and updated whenever additional data is available. For this purpose, the marginal likelihood p(y|X), which means the probability to observe the output y given the input X, is introduced. The marginal likelihood of the observations with the prior inference on the weights vector is expressed as the probability density function of the multivariate Gaussian distribution as follows: Then, the marginal likelihood will be maximized with the best combination of the hyperparameters, which yields an optimization problem. In this study, for numerical convenience, the minus log-marginal likelihood is utilized as the objective function (L): The best hyperparameters (θ best ) can be determined through optimization to minimize the objective function:

GPR with FE Analysis Results
In several previous studies, the prior mean function in GPR is set to a zero matrix, as described in Equation (4). However, the results of FE analysis can be introduced as the prior mean m(·), also called the explicit basis function, and the input-output relationship is expressed as: Compared with Equation (2), Equation (27) is seen to include the prior mean based on FE analysis, and this can help to overcome the limitation on predictions of the vertical deflection of a bridge under construction, when sufficient data is not available.
For computational simplicity, the residual differences between the output values and the explicit basis function values are considered as new outputs as follows, where y is the set of new training outputs, and f '(·) is the new regression function. Then, with the new training outputs (y ), we have an input-output relationship with a zero prior mean again, which was discussed in Section 2.1. The marginal likelihood can now be expressed as the conditional probability of the new training outputs given the training input matrix (X) as follows: Once the hyperparameters are optimized, the predictive outputs can be estimated as follows, which is equivalent to:

Approximation Algorithm for Reducing Computational Cost
One of the limitations involved in the application of the GPR is the computational cost. GPR requires the calculation of the inverse matrix of the covariance matrix, as explained in Section 2.1; hence, the cost of GPR is sensitive to the size of the input and output datasets. To address this issue, various approximation methods were developed, which use a partial dataset to derive performance that is comparable to a standard GPR, where a full dataset is applied [31,34]. Among these methods, the Subset of Data (SoD) approximation, which is simple and known to have a better speed-accuracy trade-off than other approximations [35,36], is introduced for this study.
The key idea of the SoD approximation is to select the subset that provides a reasonable hyperparameter optimization, compared with the regression, using the whole dataset. Utilizing only the M D number among the N D quantities of the full dataset (M D < N D ) naturally decreases the computational cost. The important task here is to select subset points, which facilitate an acceptable regression performance for building the posterior distribution, called the active set. Data selection algorithms, which use an active learning criterion (e.g., differential entropy [37], matching pursuit [38], or information gain [39]) have been used as more effective approaches because randomly selected subset points generally do not perform well [31,36]. This study introduces the differential entropy score suggested by Lawrence et al. (2003) as the criterion, which is maximized at the point with the largest variance [37]. Using the selection criterion, a portion of the whole dataset is selected as the active set (D A ), where X A and y A are the selected training input matrix, and output vector (i.e., active set), respectively, whose sizes are M D by N x and M D by 1. Once the active set is selected, the hyperparameters are optimized by utilizing the maximum likelihood concept described in Equations (29)-(31) with the active set. The marginal likelihood of the training output vector given the training input matrix of the active set can be described as follows, where K A is the covariance matrix of the training input matrix of the active set X A , whose size is M D by M D . The computational cost for calculating the inverse matrix in the likelihood function is much lower because of the small size of the covariance matrix. After optimizing the hyperparameters, the predictive outputs can be estimated: where K A* is the covariance matrix between the selected training (X A ) inputs and the N P number of prediction (X * ) inputs, and K A** is the covariance matrix of the prediction input matrix X * , whose sizes are M D by N P and N P by N P , respectively. When the explicit basis function is incorporated as the prior mean, the new outputs are calculated as:

Example Bridge Description
To test the proposed method, it was applied to a full-scale railway bridge under construction. The bridge consists of four PSC girders (G1, G2, G3, and G4, as shown in Figure 1a), and each girder has a trivet-shaped cross section, as shown in Figure 1b, and a width, height, and length of 2.0 m, 2.2 m, and 39.94 m, respectively. Among these four girders, G1 was selected as the target, for which actual measurements were taken, by considering various physical properties, such as vertical deflection and temperature. By following the schedule listed in Table 1, the girder experienced construction loads by the tensioning and detensioning of tendons and substructure placements. The bridge construction started with girder concrete casting and tendon tensioning, which caused an upward deflection on 18 September 2017, defined as day index 1. Tendon releasing and girder erection proceeded in By following the schedule listed in Table 1, the girder experienced construction loads by the tensioning and detensioning of tendons and substructure placements. The bridge construction started with girder concrete casting and tendon tensioning, which caused an upward deflection on 18 September 2017, defined as day index 1. Tendon releasing and girder erection proceeded in sequence on 22 September (day index 5), and 6 November (day index 50), respectively. Figure 2a shows the slab placement, which was performed from 25 January to 8 February 2018, and it caused a downward deflection. Figure 2b shows the placement of protective walls, which was performed from 14 February to 16 March 2018, and it caused additional downward deflection of the target bridge girder.

Prior prediction based on FE analysis
For the prior prediction of the time-dependent vertical deflection of the target bridge, in this study, the three-dimensional FE model, established by  is introduced [40]. The FE model was constructed using commercial FE analysis software, MIDAS/CIVIL 2017, which supports the consideration of a construction schedule and the time-dependent structural responses of creep and shrinkage [40,41]. The FE model has 212 nodes and 208 elements, and the constructed model is shown in Figure 3. In addition, the material properties of concrete and steel strands are summarized in Tables 2 and 3, respectively. The concrete for the slab, cross beam, and girder were modeled using the C30, C30, and C49 models of MIDAS/CIVIL 2017, respectively, and the SWPC7B 15.2 mm, Low Relaxation model of MIDAS/CIVIL 2017 was used for strands. The construction stages and the associated schedule in Table 1 were considered by conducting a construction stage analysis in the software, and the KCI-USD12 code [10] was utilized to estimate the time-dependent deflection of the example bridge. The FE analysis result for the vertical deflection of the bridge, during the construction stages, was introduced as the explicit basis function in the GPR, and the prior prediction was updated using actual measurement data.

Prior prediction based on FE analysis
For the prior prediction of the time-dependent vertical deflection of the target bridge, in this study, the three-dimensional FE model, established by  is introduced [40]. The FE model was constructed using commercial FE analysis software, MIDAS/CIVIL 2017, which supports the consideration of a construction schedule and the time-dependent structural responses of creep and shrinkage [40,41]. The FE model has 212 nodes and 208 elements, and the constructed model is shown in Figure 3. In addition, the material properties of concrete and steel strands are summarized in Tables 2  and 3, respectively. The concrete for the slab, cross beam, and girder were modeled using the C30, C30, and C49 models of MIDAS/CIVIL 2017, respectively, and the SWPC7B 15.2 mm, Low Relaxation model of MIDAS/CIVIL 2017 was used for strands. The construction stages and the associated schedule in Table 1 were considered by conducting a construction stage analysis in the software, and the KCI-USD12 code [10] was utilized to estimate the time-dependent deflection of the example bridge. The FE analysis result for the vertical deflection of the bridge, during the construction stages, was introduced as the explicit basis function in the GPR, and the prior prediction was updated using actual measurement data.

Measurement Data
The vertical deflection was measured using the computer vision system proposed by Lee et al. [42]. A vision system with the hardware configurations shown in Figure 4a was deployed to monitor the displacement of the main target. In this dual-camera system, the displacement of the main target can be measured, regardless of camera movement, by placing a sub-target displacement at a stationary reference point. As shown in Figure 4b, the main target was deployed at the mid-span of the girder at which displacement was monitored. With the specifications of the camera and the lens listed in Table 4, the employed computer vision system can measure displacement with a resolution of 0.2 mm at a distance of 20 m.

Measurement Data
The vertical deflection was measured using the computer vision system proposed by Lee et al. [42]. A vision system with the hardware configurations shown in Figure 4a was deployed to monitor the displacement of the main target. In this dual-camera system, the displacement of the main target can be measured, regardless of camera movement, by placing a sub-target displacement at a stationary reference point. As shown in Figure 4b, the main target was deployed at the mid-span of the girder at which displacement was monitored. With the specifications of the camera and the lens listed in Table 4, the employed computer vision system can measure displacement with a resolution of 0.2 mm at a distance of 20 m. irder at which displacement was monitored. With the specifications of the camera and the lens in Table 4, the employed computer vision system can measure displacement with a resolution mm at a distance of 20 m.  A total of 34,468 measurements were obtained for about nine months, beginning on 4 January 2018. Figure 5a shows the measurement data and the FE analysis results for the vertical deflection. Due to a malfunction in the dual-camera system, sensor data were not obtained from 10 May (day index 235) to 25 June (day index 281) 2018. The FE analysis results (represented by the green line) considers the construction schedule and the associated physical phenomena, and thus, it follows the trend of the measurements by the dual-camera system represented by black dots. However, the FE analysis can provide only a long-term deterministic prediction, without daily fluctuations, as shown in Figure 5b, which is a magnification of a portion of the plot in Figure 5a. Moreover, it is observed that the FE analysis results have a bias. That is, the FE analysis result cannot exactly predict the vertical deflection alone due to its model error.
(a) (b) Figure 5. Plots of (a) sensor data of vertical deflection versus the prior prediction from the FE analysis; (b) a magnified view of a portion of the plot in Figure 5(a).  A total of 34,468 measurements were obtained for about nine months, beginning on 4 January 2018. Figure 5a shows the measurement data and the FE analysis results for the vertical deflection. Due to a malfunction in the dual-camera system, sensor data were not obtained from 10 May (day index 235) to 25 June (day index 281) 2018. The FE analysis results (represented by the green line) considers the construction schedule and the associated physical phenomena, and thus, it follows the trend of the measurements by the dual-camera system represented by black dots. However, the FE analysis can provide only a long-term deterministic prediction, without daily fluctuations, as shown in Figure 5b, which is a magnification of a portion of the plot in Figure 5a. Moreover, it is observed that the FE analysis results have a bias. That is, the FE analysis result cannot exactly predict the vertical deflection alone due to its model error.  A total of 34,468 measurements were obtained for about nine months, beginning on 4 January 2018. Figure 5a shows the measurement data and the FE analysis results for the vertical deflection. Due to a malfunction in the dual-camera system, sensor data were not obtained from 10 May (day index 235) to 25 June (day index 281) 2018. The FE analysis results (represented by the green line) considers the construction schedule and the associated physical phenomena, and thus, it follows the trend of the measurements by the dual-camera system represented by black dots. However, the FE analysis can provide only a long-term deterministic prediction, without daily fluctuations, as shown in Figure 5b, which is a magnification of a portion of the plot in Figure 5a. Moreover, it is observed that the FE analysis results have a bias. That is, the FE analysis result cannot exactly predict the vertical deflection alone due to its model error.

Probabilistic Prediction Results
(a) (b) Figure 5. Plots of (a) sensor data of vertical deflection versus the prior prediction from the FE analysis; (b) a magnified view of a portion of the plot in Figure 5(a).

Probabilistic Prediction Results
The proposed method is applied to build probabilistic prediction models, and the flow chart, shown in Figure 6, summarizes the procedure. First, data are processed with the measurement data according to the FE model. Second, the SoD approximation is applied to reduce the computational cost. In this application example, a convergence test of the SoD approximation was conducted for the full set of 1690 data, and Figure 7 shows the results. It is observed from the figure that the root mean square errors (RMSEs) of predictive mean and variance decrease with the increasing amount of considered data. Based on these results, the minimum number of SoD was selected to be 400. Next, kernel hyperparameters were determined through an optimization procedure using SoD, and the predictive mean and variance were calculated. The proposed method is applied to build probabilistic prediction models, and the flow chart, shown in Figure 6, summarizes the procedure. First, data are processed with the measurement data according to the FE model. Second, the SoD approximation is applied to reduce the computational cost. In this application example, a convergence test of the SoD approximation was conducted for the full set of 1690 data, and Figure 7 shows the results. It is observed from the figure that the root mean square errors (RMSEs) of predictive mean and variance decrease with the increasing amount of considered data. Based on these results, the minimum number of SoD was selected to be 400. Next, kernel hyperparameters were determined through an optimization procedure using SoD, and the predictive mean and variance were calculated.

Root-Mean Square
Error (RMSE) The proposed method is applied to build probabilistic prediction models, and the flow chart, shown in Figure 6, summarizes the procedure. First, data are processed with the measurement data according to the FE model. Second, the SoD approximation is applied to reduce the computational cost. In this application example, a convergence test of the SoD approximation was conducted for the full set of 1690 data, and Figure 7 shows the results. It is observed from the figure that the root mean square errors (RMSEs) of predictive mean and variance decrease with the increasing amount of considered data. Based on these results, the minimum number of SoD was selected to be 400. Next, kernel hyperparameters were determined through an optimization procedure using SoD, and the predictive mean and variance were calculated.

Root-Mean Square
Error (RMSE) Figure 7. Convergence test result for the training data number using the SoD approximation.
The measurement data from day index 109 to 140 is used, and the prediction results with, and without, the use of the FE analysis results are presented in Figure 8a,b respectively. It is seen in both of the figures that, until day 145, the actual measurement data on the vertical deflection are primarily located within the 99% prediction interval (PI). In Figure 8a, however, the measurement data are outside the prediction interval on day 145, when the slab was placed. This is because the slab placement produced a large downward deflection on the target girder, and it could not be considered without the FE analysis. As a result, the predictive mean and standard deviation of the vertical deflection are estimated to be 20.849, and 2.099 mm, respectively. However, the actual measurement on day 145 was 12.64 mm. Using the normal distribution assumption on the vertical deflection in GPR, the probability of observing at least the level of deviation from the mean is calculated to be 2×Φ(−|12.64−20.849|/2.099) = 9.19 × 10 −5 , where Φ(·) is the standard normal cumulative distribution function (CDF), which indicates a very rare event. On the other hand, when the FE analysis result is used for the probability prediction, the predictive mean and standard deviation are estimated to be 12.639 and 2.364 mm for day 145, and the probability of observing at least the level of deviation from the mean is calculated to be 2×Φ(−|12.64−12.639|/2.364) = 0.9996, as shown in Figure 8b, which indicates a very probable event. These probabilities can be used for effective bridge management. placement produced a large downward deflection on the target girder, and it could not be considered without the FE analysis. As a result, the predictive mean and standard deviation of the vertical deflection are estimated to be 20.849, and 2.099 mm, respectively. However, the actual measurement on day 145 was 12.64 mm. Using the normal distribution assumption on the vertical deflection in GPR, the probability of observing at least the level of deviation from the mean is calculated to be 2×Φ(−|12.64−20.849|/2.099) = 9.19 × 10 −5 , where Φ(·) is the standard normal cumulative distribution function (CDF), which indicates a very rare event. On the other hand, when the FE analysis result is used for the probability prediction, the predictive mean and standard deviation are estimated to be 12.639 and 2.364 mm for day 145, and the probability of observing at least the level of deviation from the mean is calculated to be 2×Φ(−|12.64−12.639|/2.364) = 0.9996, as shown in Figure 8b, which indicates a very probable event. These probabilities can be used for effective bridge management. Similarly, predictions are made using measurement data for longer periods, as described in Figure 9, by employing both, the measurement data and the prior prediction from the FE analysis. According to a previous study [2], measurement data, obtained for approximately three months, were necessary for an accurate prediction on bridge deflection when considering environmental effects. If sufficient training data are not available, environmental effects can be considered in the prediction through their introduction into the FE analysis. However, this is not within the scope of this study. The prediction model in Figure 10a is constructed using the measurement data between the day indices 109 and 180, and the model in Figure 10b utilized the data between the day indices 109 and 287. First, Figure 10a, which describes the deflection after the placement of the protective walls shows that the vertical deflection after the protective wall placement is within the expected range. Similarly, Figure 10b shows that the deflections are within the acceptable ranges, which means that the construction has proceeded as expected. Finally, the long-term deflection prediction model was constructed using the entire measurements from day indices 109 to 330, as shown in Figure 10c. Similarly, predictions are made using measurement data for longer periods, as described in Figure 9, by employing both, the measurement data and the prior prediction from the FE analysis. According to a previous study [2], measurement data, obtained for approximately three months, were necessary for an accurate prediction on bridge deflection when considering environmental effects. If sufficient training data are not available, environmental effects can be considered in the prediction through their introduction into the FE analysis. However, this is not within the scope of this study. The prediction model in Figure 10a is constructed using the measurement data between the day indices 109 and 180, and the model in Figure 10b utilized the data between the day indices 109 and 287. First, Figure 10a, which describes the deflection after the placement of the protective walls shows that the vertical deflection after the protective wall placement is within the expected range. Similarly, Figure 10b shows that the deflections are within the acceptable ranges, which means that the construction has proceeded as expected. Finally, the long-term deflection prediction model was constructed using the entire measurements from day indices 109 to 330, as shown in Figure 10c.   As described above, the proposed method can be utilized to calculate the observance probability of an actual measurement. If the observance probability is estimated to be low, it may indicate that there is a serious problem with the bridge, and the bridge manager may consider shutting down the bridge immediately. On the other hand, if the observance probability is estimated to be high, despite the large deflection, the bridge manager can be confident that there is no problem with the bridge. In this manner, the proposed method is expected to support the efficient monitoring and management of bridges.

Conclusions
This study suggested a new nonparametric Bayesian method, which employs both an FE model and actual measurement data to construct the probabilistic prediction model of bridge deflection. To overcome the limitations of an imperfect FE model and a shortage of data, the GPR was introduced and modified to consider both, the FE analysis results and actual measurement data. In addition, the probabilistic prediction model can be updated whenever additional measurement data is available. The proposed method was applied to an actual railway bridge in the Republic of Korea, which is a PSC bridge under construction, and the probabilistic prediction results agreed with the measurements. In addition, the proposed method calculates the observance probability for a given measurement, which can be useful for efficient monitoring and management of bridges.

Conflicts of Interest:
The authors declare no conflict of interest.