Optimal Sensor Placement in Buildings: Stationary Excitation

1. Introduction

The California Strong Motion Instrumentation Program (CSMIP) in the California Geological Survey (CGS) was established in 1972, after the devastating 1971 San Fernando earthquake to obtain vital earthquake data and improve understanding of the engineering and scientific communities about earthquakes and also how civil structures respond to ground shaking [1,2,3,4,5]. Currently, the program has more than 1,390 active stations throughout California including 946 ground-response stations, 273 buildings, 26 dams, and 82 bridges with the total number of sensors exceeding 10,200. While post-earthquake damage assessment was one of its original objectives, Structural Health Monitoring (SHM) was in its infancy during CSMIP's initiation. Consequently, buildings, especially mid- to high-rise structures, were—and, to a large extent, still are —equipped with a limited number of sensors due to the sensor installation and maintenance cost (see Figure 1). Although more extensive instrumentation is recommended in some guidelines, such as the Los Angeles Tall Buildings Guidelines [6], it is not yet feasible to install sensors on every floor of a building due to instrumentation and maintenance costs.

Recent advancements over the past two decades in both software and hardware have made vibration-based post-event (e.g., earthquake) damage assessment and continuous structural health monitoring (SHM) not only feasible but essential for advancing a resilient built environment [7,8,9,10]. These capabilities were previously unattainable through traditional, labor-intensive visual in-situ inspections. Among the most common approaches to post-event assessment are data-driven modal identification methods [11,12,13,14,15,16] and model-based system identification techniques [17,18,19,20,21,22], which utilize data from instrumented structures to evaluate structural health and performance either after an earthquake or during regular operation. However, as previously noted, most instrumented structures are sparsely instrumented due to the high costs associated with installation and maintenance. This limited instrumentation introduces substantial uncertainty into response estimation [23,24,25].

Deterministic interpolation techniques, such as cubic splines [26], have been used to estimate structural responses at non-instrumented floors [27,28,29]. To address the inability of these methods to quantify uncertainty systematically, we recently introduced a hybrid deterministic-probabilistic approach [30]. This method combines a simple beam model—easily calibrated by adjusting a few parameters [31,32,33,34,35,36,37,38]—with a Gaussian Process Regression (GPR) model [39,40], a non-parametric probabilistic technique widely used across various fields, including earthquake engineering (e.g., [41]). Through extensive verification and real-world case studies, this approach has been shown to yield highly accurate mean estimations while simultaneously quantifying associated uncertainties—an essential feature when synthetic response data are intended for probabilistic post-event damage assessment, such as FEMA P-58 [42] methodologies [43].

Following the development of this hybrid approach [30], we initiated a new project aimed at minimizing estimation and interpolation uncertainty by optimizing sensor placement along the height of building structures. While this approach may not enhance existing instrumented structures, CGS has recently launched an extensive initiative to either fully re-instrument or add additional sensors to more than 70 buildings—making the findings of this sensor placement study especially relevant. Furthermore, the number of CSMIP-instrumented buildings continues to grow, as California’s seismic safety code requires seismic monitoring systems in tall buildings [6].

Optimal sensor placement has long been a subject of research in fields where data collection is essential [44,45,46,47,48,49]. In earthquake and structural engineering, interest in this topic dates back to the early years of structural vibration monitoring (see, e.g., [50,51]) and was revitalized around the year 2000 by pioneering work led by Beck and his former students [52,53,54,55,56,57]. With the exception of Limongelli’s study [58], which used a deterministic approach to minimize the difference between the response reconstructed via cubic spline interpolation and the true response for a series of typical building models, information theory forms the core of all these optimal sensor placement strategies, with differences among studies arising from the specific information metric used. Examples include the Fisher Information Matrix [50], information entropy [52], joint entropy [59], value of information [60,61], Kullback–Leibler divergence [62,63], and mutual information [56].

In addition to the information metric, the objective of the sensor placement significantly influences the results. Common objectives in the literature include finite element model updating, parameter estimation, model selection, modal identification, damage detection, and structural health monitoring. Although response reconstruction has been considered among these objectives (e.g., [57]), it is often treated as secondary to model updating. However, CSMIP’s primary mandate is to provide structural response data following earthquakes. As such, response measurement and reconstruction are central objectives.

Building on the work of Papadimitriou et al. [54] and Gorodesky and Marzouk [64], we propose an optimal sensor placement strategy for building structures that minimizes uncertainty in response prediction at non-instrumented floors using Gaussian Process Regression (GPR). The methodology is based on two key assumptions: (1) the excitation is stationary, and (2) the building is accurately modeled as a shear-flexural coupled beam. These assumptions allow us to apply random vibration theory in conjunction with beam theory to analytically construct the covariance kernel function used in the GPR model.

In the next section, we first provide a brief review of GPR. Then, using the beam model, we derive the analytical kernel function for a clamped shear-flexural beam, representative of tall building structures. Finally, we introduce the objective function used to determine the optimal sensor locations along the height of the building. A series of case studies is conducted to identify and verify optimal sensor placements for buildings with varying numbers of sensors. These results are compiled into a practical reference table. We also validate our proposed sensor placement strategy using earthquake data from an actual instrumented building. The major findings of this study are summarized in the Conclusion section.

2. The Proposed Approach

2.1. Gaussian Process Regression

Gaussian Process Regression (GPR) is a Bayesian approach for function approximation. Instead of providing a single functional value $f\left(x\right)$, GPR yields posterior probability distribution over all possible functions [40]. By definition, a GP consists of a set of random variables, any finite subset of which follows a joint Gaussian distribution. Consequently, its prior can be characterized by a prior mean function $m_0\left(x\right)$ and a covariance kernel $k\left(x,x^{\text{'}};\mathit{\vartheta }\right)$, as follows

$$f\left(x\right)~GP\left(m_0\left(x\right),k\left(x,x^{\text{'}};\mathit{\vartheta }\right)\right),$$

(1)

with

$$m_0\left(x\right)=\mathbb{E}\left[f\left(x\right)\right],$$

(2)

$$k\left(x,x^{\text{'}};\mathit{\vartheta }\right)=\mathbb{E}\left[\left(f\left(x\right)-m_0\left(x\right)\right)\left(f\left(x^{\text{'}}\right)-m_0\left(x^{\text{'}}\right)\right)\right],$$

(3)

where $\mathbb{E}\left[.\right]$ denotes the expected value, and $\mathit{\vartheta }$ is a vector of hyperparameters that define the covariance kernel. In this study, we assume that both the parameter $x$ and function value $f\left(x\right)$ are scalar. Additionally, while not necessary, the prior mean is often assumed to be zero, as any known nonzero trend can be subtracted from the data.

Now, assume that noisy function values are available at $n_d$ parameter values, given by $\overline{\mathit{x}}={\left[\begin{array}{cccc}{\overline{x}}_1& {\overline{x}}_2& \dots & {\overline{x}}_{n_d}\end{array}\right]}^T$. That is, $\overline{\mathit{y}}=\mathit{f}\left(\overline{\mathit{x}}\right)+\mathit{ϵ}$ where $\mathit{f}\left(\overline{\mathit{x}}\right)={\left[\begin{array}{cccc}f\left({\overline{x}}_1\right)& f\left({\overline{x}}_2\right)& \dots & f\left({\overline{x}}_{n_d}\right)\end{array}\right]}^T$ is the vector of true function values, $\overline{\mathit{y}}={\left[\begin{array}{cccc}{\overline{y}}_1& {\overline{y}}_2& \dots & {\overline{y}}_{n_d}\end{array}\right]}^T$ is the vector of noisy measurement, and $\mathit{ϵ}={\left[\begin{array}{cccc}{ϵ}_1& {ϵ}_2& \dots & {ϵ}_{n_d}\end{array}\right]}^T$ is the vector of measurement noise. Based on the aforementioned definition of a GP, function value at any arbitrary input value of $\widehat{x}$ along with the measurement data construct a joint Gaussian distribution. Thus,

$$\left[\begin{array}{c}\overline{\mathit{y}}\\ f\left(\widehat{x}\right)\end{array}\right]~N\left(\left[\begin{array}{c}0\\ 0\end{array}\right],\begin{array}{cc}\mathbf{C}\left(\overline{\mathit{x}},{\overline{\mathit{x}}}^{\text{'}};\vartheta \right)+{\sigma }_0^2\mathbf{I}& \mathit{k}\left(\overline{\mathit{x}},\widehat{x};\mathit{\vartheta }\right)\\ {\mathit{k}\left(\overline{\mathit{x}},\widehat{x};\mathit{\vartheta }\right)}^T& k\left(\widehat{x},\widehat{x};\mathit{\vartheta }\right)\end{array}\right),$$

(4)

where the matrix $\mathbf{C}\left(\overline{\mathit{x}},{\overline{\mathit{x}}}^{\text{'}};\mathit{\vartheta }\right)$ is the covariance matrix which is calculated using the covariance kernel as

$$\mathbf{C}\left(\overline{\mathit{x}},{\overline{\mathit{x}}}^{\text{'}};\mathit{\vartheta }\right)=\left[\begin{array}{cccc}k\left({\overline{x}}_1,{\overline{x}}_1;\mathit{\vartheta }\right)& k\left({\overline{x}}_1,{\overline{x}}_2;\mathit{\vartheta }\right)& \dots & k\left({\overline{x}}_1,{\overline{x}}_{n_d};\mathit{\vartheta }\right)\\ k\left({\overline{x}}_2,{\overline{x}}_1;\mathit{\vartheta }\right)& k\left({\overline{x}}_2,{\overline{x}}_2;\mathit{\vartheta }\right)& ⋮& k\left({\overline{x}}_2,{\overline{x}}_{n_d};\mathit{\vartheta }\right)\\ ⋮& \dots & \ddots & ⋮\\ k\left({\overline{x}}_{n_d},{\overline{x}}_1;\mathit{\vartheta }\right)& k\left({\overline{x}}_{n_d},{\overline{x}}_2;\mathit{\vartheta }\right)& \dots & k\left({\overline{x}}_{n_d},{\overline{x}}_{n_d};\mathit{\vartheta }\right)\end{array}\right],$$

(5)

and $\mathit{k}\left(\overline{\mathit{x}},\widehat{x};\mathit{\vartheta }\right)$ is a vector representing the correlation between $\mathit{f}\left(\overline{\mathit{x}}\right)$ and $f\left(\widehat{x}\right)$, that is

$$\mathit{k}\left(\mathit{x},\widehat{x};\mathit{\vartheta }\right)=\left[\begin{array}{c}k\left(x_1,\widehat{x};\mathit{\vartheta }\right)\\ k\left(x_2,\widehat{x};\mathit{\vartheta }\right)\\ ⋮\\ k\left(x_{n_d},\widehat{x};\mathit{\vartheta }\right)\end{array}\right].$$

(6)

The ${\sigma }_0^2{\mathbf{I}}_{n_d\times n_d}$ term added to the covariance matrix of the measured data is representative of the measurement noise which is assumed to be a zero-mean uncorrelated white noise with a variance of ${\sigma }_0^2$.

It is then mathematically straightforward to show that the conditional distribution of $f\left(\widehat{x}\right)$ given measurements ($\overline{\mathit{x}}$, $\overline{\mathit{y}}$) is

$$\left.f\left(\widehat{x}\right)\right|\overline{\mathit{x}},\overline{\mathit{y}}~N\left(m_{f\left(\widehat{x}\right)},{\sigma }_{f\left(\widehat{x}\right)}^2\right),$$

(7)

where

$$m_{f\left(\widehat{x}\right)}={\mathit{k}\left(\overline{\mathit{x}},\widehat{x};\mathit{\vartheta }\right)}^T{\left[\mathbf{C}\left(\overline{\mathit{x}},{\overline{\mathit{x}}}^{\text{'}};\mathit{\vartheta }\right)+{\sigma }_0^2\mathbf{I}\right]}^{-1}\overline{\mathit{y}},$$

(8)

$${\sigma }_{f\left(\widehat{x}\right)}^2=k\left(\widehat{x},\widehat{x};\mathit{\vartheta }\right)-{\mathit{k}\left(\overline{\mathit{x}},\widehat{x};\mathit{\vartheta }\right)}^T{\left[\mathbf{C}\left(\overline{\mathit{x}},{\overline{\mathit{x}}}^{\text{'}};\mathit{\vartheta }\right)+{\sigma }_0^2\mathbf{I}\right]}^{-1}\mathit{k}\left(\overline{\mathit{x}},\widehat{x};\mathit{\vartheta }\right),$$

(9)

are posterior mean and variance of the prediction at $\widehat{x}$. To apply the equations above, in addition to the data, the covariance kernel and its hyperparameters are required. The common practice is to select an appropriate covariance kernel from a wide range of available functions and then find the optimal hyperparameter values by minimizing the negative log-likelihood of the data as shown below [65]

$$NLL\left(\overline{\mathit{y}}\right)={\displaystyle \frac{1}{2}}{\overline{\mathit{y}}}^T{\left[\mathbf{C}\left(\overline{\mathit{x}},{\overline{\mathit{x}}}^{\text{'}};\mathit{\vartheta }\right)+{\sigma }_0^2\mathbf{I}\right]}^{-1}\overline{\mathit{y}}+{\displaystyle \frac{1}{2}}\log \left|\left(\mathbf{C}\left(\overline{\mathit{x}},{\overline{\mathit{x}}}^{\text{'}};\mathit{\vartheta }\right)+{\sigma }_0^2\mathbf{I}\right)\right|+{\displaystyle \frac{n_d}{2}}\log 2\pi .$$

(10)

As an example, Figure 2 compares the Squared Exponential (SE) and Matern 3/2 kernels [66], which are defined by the following equations,

$$k_{SE}\left(x,x^{\text{'}};\mathit{\vartheta }\right)={\sigma }_f^2{e}^{-{\displaystyle \frac{d^2}{2{\sigma }_l^2}}},$$

(11)

$$k_{M{\displaystyle \frac{3}{2}}}\left(x,x^{\text{'}};\mathit{\vartheta }\right)={\sigma }_f^2\left(1+{\displaystyle \frac{\sqrt{3}d}{{\sigma }_l}}\right)e^{-{\displaystyle \frac{\sqrt{3}d}{{\sigma }_l}}}.$$

(12)

where $d=\left|x-x^{\text{'}}\right|$ and $\mathit{\vartheta }={\left[{\sigma }_f^2,{\sigma }_l\right]}^T$, in which ${\sigma }_f^2$ and ${\sigma }_l$ are two hyperparameters called the signal variance and correlation length, respectively. As seen in Figure 2, assuming a signal variance of 1, a small correlation length leads to a rapid drop in correlation, indicating minimal similarity for distances greater than this length. Conversely, a large correlation length implies nearly identical function values across the domain.

For the optimal sensor placement task, data is not usually available a priori to use with Eq. (10). Therefore, we represent the structure with a simple model and use an analytical solution to determine the covariance kernel, making certain assumptions about the statistical properties of the external loads. This approach was previously used by Papadimitriou et al. [54] with a flexural beam under simply supported boundary conditions, which are not applicable to building structures. In this study, we reproduce the kernel function using the cantilever coupled beam model [31,32,33,67] under spatially deterministic (uniform) but temporally uncorrelated lateral loads.

2.2. Analytical Covariance Kernel Using Beam Model

Figure 3a shows the beam model used in this study, which combines flexural and shear beams, shown in Figure 3b and Figure 3c, respectively. Assuming a uniform mass per unit length $\rho$, shear stiffness ($GA$), and flexural stiffness ($EI$) along the beam height $L$, the governing partial differential equation of this model under lateral force $p\left(x,t\right)$ is given as [31]

$$\rho {\displaystyle \frac{{\partial }^{2}u\left(x,t\right)}{\partial t^2}}+{\displaystyle \frac{EI}{L^4}}{\displaystyle \frac{{\partial }^{4}u\left(x,t\right)}{\partial x^4}}-{\displaystyle \frac{GA}{L^2}}{\displaystyle \frac{\partial u\left(x,t\right)}{\partial x}}=p\left(x,t\right),$$

(13)

where $u\left(x,t\right)$ represents lateral displacement response of the beam at normalized height $x=h/L$ ($0\le x\le 1$). Dividing both sides of the equation above by $EI$, we can write

$${\displaystyle \frac{\rho }{EI}}{\displaystyle \frac{{\partial }^{2}u\left(x,t\right)}{\partial t^2}}+{\displaystyle \frac{1}{L^4}}{\displaystyle \frac{{\partial }^{4}u\left(x,t\right)}{\partial x^4}}-{\displaystyle \frac{{\alpha }^2}{L^4}}{\displaystyle \frac{\partial u\left(x,t\right)}{\partial x}}={\displaystyle \frac{1}{EI}}p\left(x,t\right),$$

(14)

where $\alpha =L\sqrt{{\displaystyle \frac{GA}{EI}}}$ is a dimensionless parameter that indicates the degree of relative contributions of flexural and shear behaviors. Using the method of separation of variables, where $u\left(x,t\right)={\sum }_{k=1}^{\infty }{\phi }_k\left(x\right)q_k\left(t\right)$, the equation of motion under free vibration conditions can be written as

$${\displaystyle \frac{EI}{\rho L^4}}\left({\displaystyle \frac{{\phi }_k^{\text{'}\text{'}\text{'}\text{'}}\left(x\right)}{{\phi }_k\left(x\right)}}-{\alpha }^2{\displaystyle \frac{{\phi }_k^{\text{'}\text{'}}\left(x\right)}{{\phi }_k\left(x\right)}}\right)=-{\displaystyle \frac{{\ddot{q}}_k\left(t\right)}{q_k\left(t\right)}}={\omega }_k^2,$$

(15)

resulting in two sets of separate equations of

$${\ddot{q}}_k\left(t\right)+{\omega }_k^2{q}_k\left(t\right)=0,$$

(16)

$${\phi }_k^{\text{'}\text{'}\text{'}\text{'}}\left(x\right)-{\alpha }^2{\phi }_k^{\text{'}\text{'}}\left(x\right)-{\omega }_k^2{\displaystyle \frac{\rho L^4}{EI}}{\phi }_k\left(x\right)=0$$

(17)

Solving the characteristic equation of Eq. (17) results in the following real (${\beta }_k$) and imaginary ($i{\gamma }_k$) roots as

$${\beta }_k=\sqrt{{\displaystyle \frac{1}{2}}{\alpha }^2+{\displaystyle \frac{1}{2}}\sqrt{{\alpha }^4+4{\omega }_k^2{\displaystyle \frac{\rho L^4}{EI}}}},$$

(18)

$${\gamma }_k=\sqrt{-{\displaystyle \frac{1}{2}}{\alpha }^2+{\displaystyle \frac{1}{2}}\sqrt{{\alpha }^4+4{\omega }_k^2{\displaystyle \frac{\rho L^4}{EI}}}}.$$

(19)

From which it is easy to show that ${\beta }_k^2={\alpha }^2+{\gamma }_k^2$. Therefore, the natural frequencies of the beam can be expressed versus ${\gamma }_k$ and ${\beta }_k$ as follows

$${\omega }_k^2={\displaystyle \frac{EI}{\rho L^4}}{\gamma }_k^2{{\beta }_k}^2.$$

(20)

Having the roots of the characteristic equation, the solution of Eq. (17), i.e., mode shapes, can be written as

$${\phi }_k\left(x\right)=c_1\sin \left({\gamma }_{k}x\right)+c_2\cos \left({\gamma }_{k}x\right)+c_3\sinh \left({\beta }_{k}x\right)+c_4\cosh \left({\beta }_{k}x\right)$$

(21)

Using four boundary conditions, and solving the system of homogenous equations, ${\phi }_k\left(x\right)$ can be calculated up to an arbitrary factor $d_k$ as

$${\phi }_k\left(x\right)=d_k\left[\sin \left({\gamma }_{k}x\right)-{\displaystyle \frac{{\gamma }_k}{{\beta }_k}}\sinh \left({\beta }_{k}x\right)+{\eta }_k\left[\cosh \left({\beta }_{k}x\right)-\cos \left({\gamma }_{k}x\right)\right]\right]$$

(22)

with

$${\eta }_k={\displaystyle \frac{{\gamma }_k^2\sin {\gamma }_k+{\gamma }_k{\beta }_k\sinh {\beta }_k}{{\gamma }_k^2\cos {\gamma }_k+{\beta }_k^2\cosh {\beta }_k}},$$

(23)

and an additional equation

$$2+\left(2+{\displaystyle \frac{{\alpha }^4}{{\gamma }_k^2{{\beta }_k}^2}}\right)\cos \left({\gamma }_k\right)\cosh \left({\beta }_k\right)+{\displaystyle \frac{{\alpha }^2}{{\gamma }_k{\beta }_k}}\sin \left({\gamma }_k\right)\sinh \left({\beta }_k\right)=0,$$

(24)

from which ${\gamma }_k$ that are needed to calculate the natural frequencies in Eq. (20), can be determined. To complete the solution, temporal component of the response also needs to be calculated. Instead of Eq. (16), which applies to free vibration conditions, we solve Eq. (14), having estimated the spatial part of the response. That is,

$${\displaystyle \frac{\rho }{EI}}{\sum }_{k=1}^{\infty }{\phi }_k\left(x\right){\ddot{q}}_k\left(t\right)+{\displaystyle \frac{1}{L^4}}{\sum }_{k=1}^{\infty }\left({\phi }_k^{\text{'}\text{'}\text{'}\text{'}}\left(x\right)-{\alpha }^2{\phi }_k^{\text{'}\text{'}}\left(x\right)\right)q_k\left(t\right)={\displaystyle \frac{1}{EI}}p\left(x,t\right).$$

(25)

By multiplying both sides of the equation above by ${\phi }_j\left(x\right)$ and integrating over the length of the beam, we have

$$\rho {\sum }_{k=1}^{\infty }{\ddot{q}}_k\left(t\right){\int }_0^1{\phi }_k\left(x\right){\phi }_j\left(x\right)dx+{\displaystyle \frac{EI}{L^4}}{\sum }_{k=1}^{\infty }{q}_k\left(t\right){\int }_0^1{\phi }_j\left(x\right)\left({\phi }_k^{\text{'}\text{'}\text{'}\text{'}}\left(x\right)-{\alpha }^2{\phi }_k^{\text{'}\text{'}}\left(x\right)\right)dx={\int }_0^1{\phi }_j\left(x\right)p\left(x,t\right)dx,$$

(26)

which can be simplified as

$$M_j{\ddot{q}}_j\left(t\right)+K_j{q}_j\left(t\right)=P_j\left(t\right),\mathrm{for}j=1\dots \infty$$

(27)

where

$$M_j=\rho {\int }_0^1{\phi }_j^2\left(x\right)dx,$$

(28)

$$K_j={\displaystyle \frac{{\alpha }^{2}EI}{L^4}}{\int }_0^1{\left[{\phi }_j^{\text{'}}\left(x\right)\right]}^{2}dx,$$

(29)

$$P_j\left(t\right)={\int }_0^1{\phi }_j\left(x\right)p\left(x,t\right)dx,$$

(30)

in which we used the fact that

$${\int }_0^1{\phi }_j\left(x\right){\phi }_k\left(x\right)dx=\left\{\begin{array}{c}0\mathrm{f}\mathrm{o}\mathrm{r}j\ne k\\ {\int }_0^1{\phi }_j^2\left(x\right)dx\mathrm{f}\mathrm{o}\mathrm{r}j=k\end{array}\right.,$$

(31)

$${\int }_0^1{\phi }_j\left(x\right)\left({\phi }_k^{\mathrm{\text{'}}\mathrm{\text{'}}\mathrm{\text{'}}\mathrm{\text{'}}}\left(x\right)-{\alpha }^2{\phi }_k^{\mathrm{\text{'}}\mathrm{\text{'}}}\left(x\right)\right)=\left\{\begin{array}{c}0\mathrm{f}\mathrm{o}\mathrm{r}j\ne k\\ {\int }_0^1{\left[{\phi }_j^{\mathrm{\text{'}}}\left(x\right)\right]}^{2}dx\mathrm{f}\mathrm{o}\mathrm{r}j=k\end{array}\right..$$

(32)

For implementation convenience, the analytical expression of these two integrals is provided below

$${\displaystyle \frac{1}{d_j^2}}{\int }_0^1{\phi }_j^2\left(x\right)dx={\displaystyle \frac{1}{4}}\left(2-{\displaystyle \frac{2{\gamma }_j^2}{{\beta }_j^2}}-{\displaystyle \frac{2{\eta }_j}{{\gamma }_j}}+{\displaystyle \frac{2{\gamma }_j{\eta }_j}{{\beta }_j^2}}+4{\eta }_j^2+{\displaystyle \frac{2{\eta }_j\cos 2{\gamma }_j}{{\gamma }_j}}-{\displaystyle \frac{2{\gamma }_j{\eta }_j\cosh 2{\beta }_j}{{\beta }_j^2}}-{\displaystyle \frac{8{\gamma }_j\left(1+{\eta }_j^2\right)\cosh {\beta }_j\sin {\gamma }_j}{{\beta }_j^2+{\gamma }_j^2}}+{\displaystyle \frac{\left({\eta }_j^2-1\right)\sin 2{\gamma }_j}{{\gamma }_j}}+{\displaystyle \frac{8\left({\gamma }_j^2-{\beta }_j^2{\eta }_j^2\right)\cos {\gamma }_j\sinh {\beta }_j}{{\beta }_j\left({\beta }_j^2+{\gamma }_j^2\right)}}+{\displaystyle \frac{8{\eta }_j\sin {\gamma }_j\sinh {\beta }_j}{{\beta }_j}}+{\displaystyle \frac{\left({\gamma }_j^2+{\beta }_j^2{\eta }_j^2\right)\sinh 2{\beta }_j}{{\beta }^3}}\right)$$

(33)

$${\displaystyle \frac{1}{d_j^2}}{\int }_0^1{\left[{\phi }_j^{\mathrm{\text{'}}}\left(x\right)\right]}^{2}dx={\displaystyle \frac{1}{4}}\left(-4{\gamma }_j{\eta }_j-2{\beta }_j^2{\eta }_j^2+2{\gamma }_j^2\left(2+{\eta }_j^2\right)-2{\gamma }_j{\eta }_j\cos 2{\gamma }_j-2{\gamma }_j{\eta }_j\cosh 2{\beta }_j-{\displaystyle \frac{8{\gamma }_j\left({\gamma }_j^2-{\beta }_j^2{\eta }_j^2\right)\cosh {\beta }_j\sin {\gamma }_j}{{\beta }_j^2+{\gamma }_j^2}}-{\gamma }_j\left({\eta }_j^2-1\right)\sin 2{\gamma }_j+8{\gamma }_j\cos {\gamma }_j\left({\eta }_j\cosh {\beta }_j-{\displaystyle \frac{{\beta }_j{\gamma }_j\left(1+{\eta }_j^2\right)\sinh {\beta }_j}{{\beta }_j^2+{\gamma }_j^2}}\right)+{\displaystyle \frac{\left({\gamma }_j^2+{\beta }_j^2{\eta }_j^2\right)\sinh 2{\beta }_j}{{\beta }_j}}\right).$$

(34)

Dividing Eq. (27) by $M_j$ and adding viscous damping with modal damping ratio of ${\xi }_j$ (We add damping to the equations at this point to keep the derivations simple.), we have

$${\ddot{q}}_j\left(t\right)+2{\xi }_j{\omega }_j{\dot{q}}_j\left(t\right)+{\omega }_j^2{q}_j\left(t\right)={\displaystyle \frac{P_j}{M_j}}\left(t\right),\mathrm{for}j=1\dots \infty$$

(35)

where ${\omega }_j^2={\displaystyle \frac{K_j}{M_j}}$ is the j-th undamped natural frequency of the beam, as previously introduced in Eq. (20). As mode shapes are arbitrary scaled (see Eq. (22)), we use mass-normalized mode shapes to remove $M_j$ from equation above. To do so, we define ${\varphi }_j\left(x\right)={\displaystyle \frac{{\phi }_j\left(x\right)}{\sqrt{M_j}}}$, and rewrite Eq. (35) as

$${\ddot{q}}_j\left(t\right)+2{\xi }_j{\omega }_j{\dot{q}}_j\left(t\right)+{\omega }_j^2{q}_j\left(t\right)=P_j\left(t\right),\mathrm{for}j=1\dots \infty ,$$

(36)

in which $P_j\left(t\right)$ is calculated using Eq. (30) after replacing ${\phi }_j\left(x\right)$ with ${\varphi }_j\left(x\right)$. The solution to Eq. (36) can be obtained using the Duhamel integral [68] as $q_j\left(t\right)={\int }_0^t{h}_j\left(t-\tau \right)P_j\left(t\right)d\tau$ where $h_j\left(t\right)$ is the continuous-time Impulse Response Function (IRF), defined as

$$h_j\left(t\right)={\displaystyle \frac{1}{{\omega }_{jd}}}{e}^{-{\xi }_j{\omega }_{j}t}\sin {\omega }_{jd}t,\mathrm{for}j=1\dots \infty$$

(37)

where ${\omega }_{jd}=\sqrt{1-{\xi }_j^2}$. The final solution, i.e., the lateral displacement of the beam at level $x$ and time $t$, can be calculated as

$$u\left(x,t\right)\cong {\sum }_{j=1}^m{\varphi }_j\left(x\right)q_j\left(t\right),$$

(38)

where the upper bound is restricted to include only the first $m$ significant modes.

Using Eq. (38), correlation between beam responses at two points $x_1$ and $x_2$ can be calculated as

$$k_{u,u}\left(x_1,x_2,t\right)\cong {\sum }_{k=1}^m{\sum }_{j=1}^m{\varphi }_k\left(x_1\right){\varphi }_j\left(x_2\right)k_{q_k,q_j}\left(t,t\right),$$

(39)

where $k_{q_k,q_j}\left(t,t\right)$ is the correlation between modal coordinates $q_k\left(t\right)$ and $q_j\left(t\right)$, which can be calculated as

$$k_{q_k,q_j}\left(t,t\right)={\int }_0^t{\int }_0^t{k}_{P_k,P_j}\left({\tau }_1,{\tau }_2\right)h_k\left(t-{\tau }_1\right)h_j\left(t-{\tau }_2\right)d{\tau }_{1}d{\tau }_2,$$

(40)

in which $k_{P_k,P_j}\left({\tau }_1,{\tau }_2\right)$ is the correlation between modal loads $P_k\left(t\right)$ and $P_j\left(t\right)$. Similar to Papadimitriou et al. (2005), we assume that the lateral load can be expressed in a separable form $p\left(x,t\right)=s\left(x\right)r\left(t\right)$. Here, we further assume that the load is spatially deterministic and constant along the length of the beam, i.e., $s\left(x\right)=s_0$, so the only source of randomness arises from its time variation, $r\left(t\right)$. Assuming this temporal randomness follows a zero-mean, stationary Gaussian white noise process with a variance of ${\mathsf{\sigma }}_r^2\delta \left(t_1-t_2\right)$, $k_{P_k,P_j}\left(t_1,t_2\right)$ can be calculated as

$$k_{P_k,P_j}\left(t_1,t_2\right)={\widehat{s}}_{k,j}{\mathsf{\sigma }}_r^2\delta \left(t_1-t_2\right),$$

(41)

where ${\widehat{s}}_{k,j}=s_0^2\left({\int }_0^1{\varphi }_k\left(x\right)dx\right)\left({\int }_0^1{\varphi }_j\left(x\right)dx\right)$ in which

$$\sqrt{M_j}{\int }_0^1{\varphi }_j\left(x\right)dx={\displaystyle \frac{1-\cos {\gamma }_j}{{\gamma }_j}}-{\displaystyle \frac{{\gamma }_j}{{\beta }_j^2}}\left(\cos {\beta }_j-1\right)+{\eta }_j\left[{\displaystyle \frac{\sinh {\beta }_j}{{\beta }_j}}-{\displaystyle \frac{\sin {\gamma }_j}{{\gamma }_j}}\right].$$

(42)

Therefore, when the modal responses reach stationary conditions (for large $t$), $k_{q_k,q_j}$ can be calculated as

$$k_{q_k,q_j}={\widehat{s}}_{k,j}{\displaystyle \frac{2{\nu }_{k,j}}{{\left({\omega }_{kd}^2-{\omega }_{jd}^2\right)}^2+{\nu }_{k,j}^4+2{\nu }_{k,j}^2\left({\omega }_{kd}^2+{\omega }_{jd}^2\right)}}{\mathsf{\sigma }}_r^2,$$

(43)

where ${\nu }_{k,j}={\xi }_k{\omega }_k+{\xi }_j{\omega }_j$. Thus, the correlation between beam responses at two points $x_1$ and $x_2$ under stationary conditions becomes

$$k_{u,u}\left(x_1,x_2\right)\cong {\sum }_{k=1}^m{\sum }_{j=1}^m{\varphi }_k\left(x_1\right){\varphi }_j\left(x_2\right)k_{q_k,q_j}.$$

(44)

This covariance kernel can be used in Eqs. (8) and (9) to predict the posterior mean and variance of the response at an unmeasured point $\widehat{x}$ using information from measurement points $\overline{\mathit{x}}$, paving the way for optimal sensor placement, as will be further described in the next section.

2.3. The Objective Function

As mentioned earlier, the purpose of this study is to develop a method for determining the optimal sensor locations that result in a more reliable reconstruction of the response of non-instrumented floors. Gaussian Process Regression is a powerful non-parametric approach that can be used to achieve this objective. To maintain consistency, we determine the optimal sensor locations by minimizing the uncertainty in the GPR prediction, following the approach used by other researchers (e.g., [64]). Assuming there are $n_d$ sensors to be distributed along the height of the building, with normalized height $x\in \left[\mathrm{0,1}\right]$, the optimal location of sensors are given by ${\mathit{x}}^{*}={\left[\begin{array}{cccc}{{\overline{x}}^{*}}_1& {{\overline{x}}^{*}}_2& \dots & {{\overline{x}}^{*}}_{n_d}\end{array}\right]}^T$ and is obtained by solving the following minimization problem:

$${\mathit{x}}^{*}=\underset{\overline{\mathit{x}}\in S}{\mathit{min}}\int \left(\left.{\sigma }_{f\left(\widehat{x}\right)}^2\right|\overline{\mathit{x}}\right)d\widehat{x},$$

(45)

where $\overline{\mathit{x}}$ is the vector of instrumented locations, $S$ is the collection of all feasible/allowable points for instrumentation, and ${\sigma }_{f\left(\widehat{x}\right)}^2$ is the posterior variance of the response at non-instrumented levels calculated using Eq. (9) and the covariance kernel constructed in the previous section. For practical implementation, the objective function in Eq. (45) can be expressed as a summation as follows

$$OF={\displaystyle \frac{1}{n_0}}{\sum }_{k=1}^{n_0}\left.{\sigma }_{f\left({\widehat{x}}_k\right)}^2\right|\overline{\mathit{x}},$$

(46)

where ${\widehat{x}}_k$ for $k=1\dots n_0$ represents the normalized height of the building’s floors with $n_0$ levels. The optimization should ideally be carried out in the discrete domain; however, we keep it continuous and convert the optimal solution ${\mathit{x}}^{*}$ to real-life by selecting the nearest floors.

As shown in Eq. (9), response data is not required to compute ${\sigma }_{f\left({\widehat{x}}_k\right)}^2$; therefore, the optimal sensor locations can be determined prior to sensor deployment, which is a key advantage of the proposed objective function. The only required information pertains to the characteristics of the beam model, which are necessary to compute the covariance kernel using Eq. (44). Specifically, the parameters that contribute to the objective function include $s_0$, ${\mathsf{\sigma }}_r^2$, $\alpha$, $\rho$, ${\omega }_j$ and ${\xi }_j(\mathrm{f}\mathrm{o}\mathrm{r}j=1\dots m)$. From Eq. (20), it is evident that once the first natural frequency ${\omega }_1$ is known, all other natural frequencies can be derived as

$${\omega }_k^2={\displaystyle \frac{{\gamma }_k^2\left({\alpha }^2+{\gamma }_k^2\right)}{{\gamma }_1^2\left({\alpha }^2+{\gamma }_1^2\right)}}{\omega }_1^2,$$

(47)

in which $\alpha$ is the only required parameter. Once $\alpha$ is known, the mode shapes can also be determined up to a free scaling factor. To obtain mass-normalized mode shapes, the mass density $\rho$ is needed. However, since $\rho$ introduces a constant scaling factor across all modes, it uniformly scales $k_{u,u}\left(x_1,x_2\right)$ at all locations. The same applies to the load parameters $s_0$ and ${\mathsf{\sigma }}_r^2$ (see Eq. (43)). Therefore, we can assume that ${\widehat{s}}_{k,j}\propto \left({\int }_0^1{\varphi }_k\left(x\right)dx\right)\left({\int }_0^1{\varphi }_j\left(x\right)dx\right)$. As a result, the only remaining controlling parameters are ${\omega }_1$, $\alpha$ and ${\xi }_j(\mathrm{f}\mathrm{o}\mathrm{r}j=1\dots m)$. It is reasonable to assume that modal damping ratios are of similar order across all modes for the purpose of optimal sensor placement, particularly when no structural response data is available. Based on this assumption, we can further simplify $k_{q_k,q_j}$ as

$$k_{q_k,q_j}\propto {\widehat{s}}_{k,j}{\displaystyle \frac{1}{{\omega }_1^3}}{\displaystyle \frac{2\xi \left(b_k+b_j\right)}{\left(1-{\xi }^2\right){\left(b_k^2-b_j^2\right)}^2+{\xi }^4{\left(b_k+b_j\right)}^4+2{\xi }^2\left(1-{\xi }^2\right){\left(b_k+b_j\right)}^2\left(b_k^2+b_j^2\right)}},$$

(48)

where $b_k={\gamma }_k^2\left({\alpha }^2+{\gamma }_k^2\right)/{\gamma }_1^2\left({\alpha }^2+{\gamma }_1^2\right)$. This expression indicates that ${\omega }_1$ affects all auto- and cross-modal correlations in a similar manner and can therefore be omitted when computing from the objective function of Eq. (46). Furthermore, assuming small damping (approximately 5%), Eq. (48) can be simplified as

$$k_{q_k,q_j}\propto {\widehat{s}}_{k,j}{\displaystyle \frac{2\xi }{\left(b_k+b_j\right)\left[\left(1+2{\xi }^2\right)\left(b_k^2+b_j^2\right)-2b_k{b}_j\right]}}.$$

(49)

When $k=j$, $k_{q_k,q_k}\propto {\widehat{s}}_{k,j}{\displaystyle \frac{1}{4b_k^3\xi }}$, indicating that changes in the damping ratio affect the modal coordinate autocorrelations in a uniform manner. Therefore, we can write: $k_{q_k,q_k}\propto {\widehat{s}}_{k,k}{\displaystyle \frac{1}{4b_k^3}}$. For $k\ne j$, we can neglect the $2{\xi }^2$ term in comparison to 1 and express: $k_{q_k,q_j}\propto {\widehat{s}}_{k,j}{\displaystyle \frac{2\xi }{\left(b_k+b_j\right){\left(b_k-b_j\right)}^2}}$. Here, damping appears again as a scaling factor affecting all cross-correlations similarly, though with a different order than in the autocorrelation case. However, since modal cross-correlations are generally much smaller than autocorrelations, the damping ratio can also be excluded from the list of controlling parameters. Thus, the only remaining parameter is α, which characterizes the relative contribution of flexural and shear behavior and can be estimated based on the building’s lateral force-resisting system. For example [32]:

• In shear wall and braced frame buildings, $\alpha$ typically ranges between 0 and 1.5.
• In dual structural systems—such as a combination of moment-resisting frames with shear walls or braced frames— $\alpha$ usually falls between 1.5 and 5.
• In moment-resisting frame buildings, $\alpha$ typically ranges from 5 and 20.

Figure 4 shows the first five mode shapes of the coupled beam model with $\alpha$ ranging from 0 to 30. The blue curve represents pure flexural behavior, while the red curve represents pure shear behavior. As shown, the nodes in each mode move toward the base of the beam as the behavior transitions from flexural to shear dominance. Additionally, the differences between mode shapes diminish for higher modes.

3. Case Studies

3.1. Single Sensor

First, we present the optimal sensor placement when only one sensor is deployed. While more than one sensor is usually used in real applications, this single-sensor scenario helps to understand and verify the details of the methodology. Figure 5 shows the results for a case with $\alpha =10$. In Figure 5a, the response variance at each point along the height of the beam (${\sigma }_{f\left(\widehat{x}\right)}^2$ in Eq. (9)) is shown versus the location of the sensor. In this heatmap plot, all responses are normalized to the prior variance of the roof, $k\left(\mathrm{1,1}\right)$, to eliminate the effects of load intensity on the absolute variance values for comparison purposes. We will refer to this scaling scheme as uniform scaling in this paper. As can be seen, the variance is close to zero along the diagonal because the sensor is placed there, and it increases as the response point moves farther from the measurement points—especially when the response point is close to the roof while the sensor is placed near the base. The objective function defined in Eq. (46) is computed by averaging along the vertical direction of the heatmap plot and is shown in the top subplot. As seen, the minimum total uncertainty is obtained when the single sensor is placed at a normalized height of 0.71.

Figure 6b presents similar plots, but with each response variance scaled by its prior variance (i.e., ${\sigma }_{f\left(\widehat{x}\right)}^2/k\left(\widehat{x},\widehat{x}\right)$). In this paper, we will refer to this scaling scheme as nonuniform scaling. With this normalization, the response uncertainty changes uniformly around the sensor location. Consequently, the total uncertainty reaches its minimum when the sensor is positioned at the mid-height of the building, as shown in the top subplot of Figure 6b.

To further evaluate the performance of this optimal sensor placement strategy, we subjected the beam model with $\alpha =10$ and a 5% modal damping ratio to 500 randomly generated white noise realizations, ensuring a sufficient duration to achieve a stationary response. Assuming sensor single-sensor measurements, we then predicted the mean responses along the height of the beam for each case using Eq. (8). Finally, we calculated the Mean Squared Error (MSE) of the predictions for two scaling scenarios as follows

$$MSE\left(\overline{x}\right)={\displaystyle \frac{1}{n_0}}{\sum }_{k=1}^{n_0}{\displaystyle \frac{{\sum }_{n=1}^{n=N}{\left[u_m\left({\widehat{x}}_k,n\right)-u({\widehat{x}}_k,n)\right]}^2}{{\sum }_{n=1}^{n=N}{u(1,n)}^2}},$$

(50)

$$MSE\left(\overline{x}\right)={\displaystyle \frac{1}{n_0}}{\sum }_{k=1}^{n_0}{\displaystyle \frac{{\sum }_{n=1}^{n=N}{\left[u_m\left({\widehat{x}}_k,n\right)-u({\widehat{x}}_k,n)\right]}^2}{{\sum }_{n=1}^{n=N}{u({\widehat{x}}_k,n)}^2}},$$

(51)

where $u({\widehat{x}}_k,n)$ and $u_m\left({\widehat{x}}_k,n\right)$ are, respectively, the true and mean predicted response at the normalized level ${\widehat{x}}_k$ and time instant $n$. $N$ is the length of the signals, and $n_0$ is the total number of discrete floors, as introduced earlier.

Figure 6 shows the Probability Distribution Function (PDF) of the MSE versus the location of a single sensor along the height of the beam for the two scaling scenarios. As expected, the mean MSE is minimized when the sensor is placed at a height of 0.71 when the roof variance is used as the scaling factor. In contrast, the mid-height of the beam is the optimal sensor location when each floor’s MSE is scaled by its response variance. Additionally, both the mean MSE and its uncertainty are minimized when the sensor is at the optimal location, as seen in both plots of Figure 6. Finally, we repeat the optimal sensor placement analysis for the entire range of the parameter $\alpha$ (i.e., 0 to 30). The total uncertainty for all of these cases and both scaling scenarios is shown in Figure 7. Similar to the mode shape plots (Figure 4), cases with pure flexural and pure shear behaviors are highlighted in blue and red, respectively. As observed, when the variances are not scaled (or scaled using a constant, such as the roof’s variance), the optimal sensor location shifts from higher to lower levels as the beam behavior transitions from flexural to shear-dominated. However, under the second scaling scheme, the optimal location remains nearly unchanged and even shifts slightly upward. While changes in the optimal location versus $\alpha$ are generally negligible in this case, the total uncertainty obtained from the two most extreme cases (i.e., pure shear and pure flexural) does not perfectly represent the lower and upper envelopes of all cases. Figure 8 presents the optimal sensor location versus $\alpha$ for both scaling scenarios. As mentioned earlier and shown in this figure, when nonuniform scaling is employed, the optimal location remains nearly unchanged. However, with uniform scaling, the optimal location shifts to lower elevations as the shear contribution increases.

3.2. Multiple Sensors

When the number of sensors is greater than one, it may not be possible to graphically present some of the results shown in the previous section. Figure 9 illustrates the total uncertainty (objective function) of a beam with $\alpha =10$ versus a normalized index representing the sensor layout. For each sensor count, we calculated the objective function for all possible sensor arrangements within the available $n_0$ locations, then sorted the resulting objective values from minimum to maximum. Because the number of possible combinations varies with the number of sensors, we normalized the sorted layout indices (ranging from 1 to the total number of combinations) by the maximum number of combinations to enable plotting all cases on the same figure.

The single-sensor scenario in Figure 9 corresponds to the top plots of Figure 5, now simply sorted from minimum to maximum. As seen in Figure 9, adding a second sensor significantly reduces the objective function, regardless of the sensor locations. This reduction can reach up to two orders of magnitude if the sensors are optimally placed. The figure demonstrates that, under both scaling schemes, using two sensors instead of one results in lower total uncertainty—even if the sensors are randomly located in the building. This trend of decreasing total uncertainty continues with the addition of more sensors, as shown in Figure 9. The number of sensors is limited to four in this figure because only five modes are considered in the analysis; thus, using five sensors would eliminate the uncertainty entirely. This figure can also be used to determine both the number and locations of sensors required to achieve a specified level of total uncertainty.

Figure 10 shows the optimal sensor locations corresponding to the minimum objective function. In this figure, blue circles indicate the optimal locations under the uniform scaling scenario, while red circles represent the optimal locations under the nonuniform scaling scenario. As seen, and similar to the single-sensor case, scaling the response variances by their corresponding prior variances results in a downward shift in the optimal sensor locations.

Similar to the single-sensor scenario, we now solve the problem for a range of α values, varying from zero (pure flexural beam) to 30 (pure shear beam). Figure 11 shows the minimum total uncertainty for different numbers of sensors under both uniform and nonuniform scaling scenarios. As seen in the figure, for both scaling cases, once $\alpha$ exceeds a certain threshold (~15 for uniform scaling and ~10 for nonuniform scaling), the minimum total uncertainty becomes almost insensitive to further changes in α. In other words, the minimum total uncertainty is more sensitive to the balance between shear and flexural contributions when the beam’s behavior is dominated by flexural effects—typical of buildings with concrete shear walls.

Figure 12 shows the optimal sensor heights for various numbers of sensors under both uniform and nonuniform scaling scenarios. As observed in the single-sensor case, when the building behavior shifts toward shear behavior, the optimal sensor locations move toward the lower levels. Additionally, and consistent with the observations in Figure 10 for a fixed value of $\alpha$, nonuniform scaling further shifts the optimal sensor locations downward and arranges them in a symmetric pattern with respect to the building’s mid-height.

To make the results of this study applicable to engineering practice, we recommend optimal sensor locations for three common structural systems, based on the mean α values suggested by [32], under both scaling scenarios. In the absence of detailed information about a specific building, if the structure can be reasonably approximated by a beam model, this table can guide sensor placement along the building height to enable estimation of non-instrumented floor response using Gaussian Process Regression with minimal total uncertainty.

4. Validation

In this section, data from an instrumented building is used to validate the optimal sensor placement method described and tested in the previous sections. A complete validation ideally requires a tall building instrumented on every floor, with data that align with the load assumptions made during the method’s derivation. Unfortunately, to the best of the authors’ knowledge, no such building exists in the CSMIP database. One of the best available cases for validation is the One Rincon Hill (ORH) building in San Francisco (see Figure 13). This 62-story building features a concrete core shear wall supplemented by tall outrigger columns. The core is connected to the outrigger columns using steel buckling-restrained braces. Additionally, the building is equipped with rooftop water tanks that function as liquid mass dampers. The floors consist of post-tensioned concrete flat slabs, and the building is supported by a 12-ft-thick mat foundation resting on bedrock. The ORH building was instrumented in 2012 through a joint effort by the California Geological Survey (CGS) and the United States Geological Survey (USGS), using 72 accelerometers across 26 levels. Figure 13 illustrates the existing instrumentation, with half of the sensors managed by CGS and the other half by USGS.

While some publicly available stationary ambient data exist, we use one of the weak earthquake records available through CESMD (Center for Engineering Strong Motion Data) (strongmotioncenter.org, last updated: 5/15/2025) for validation to ensure that the applied lateral inertia forces are spatially deterministic and uniform. The 2014 M6 South Napa earthquake, which occurred approximately 50 km from the building, produced a Peak Ground Acceleration (PGA) of 0.005g. This level of ground motion satisfies the assumption of stationarity and temporal uncorrelation used in the formulation. Figure 14 shows the recorded displacement signals (NS direction) used to identify optimal sensor locations for this building under different scaling and sensor quantity scenarios. Since the building's floor plan changes significantly at Level 7, we consider that level as the base of the structure for calculating floor heights and selecting data channels. This results in a total of 23 instrumented floors used in the analysis.

To verify the findings reported in Table 1, we predict the mean response of the non-instrumented floors under each scaling and sensor quantity scenario using Eq. (8), and then compute the errors using Eqs. (50) and (51). The minimum error is expected when the sensor is placed at the optimal location. Figure 15 shows the variation of the mean-square error for the single-sensor scenario. This figure closely resembles Figure 6, although it is based on real data from a building that does not perfectly behave as a shear-flexural beam (an α value of 1.5 is used in the beam model calculations). Also, the loading is not ideal stationary white noise. As seen in the figure, the optimal sensor locations are approximately 0.83 and 0.58 times the building height for the uniform and nonuniform scaling scenarios, respectively, consistent with the results suggested in Table 1.

Figure 16 displays the identified optimal sensor locations for various numbers of sensors under the two scaling scenarios. As shown in the figure, the identified sensor locations closely match those suggested in Table 1 across different sensor quantities.

5. Conclusions

We developed an optimal sensor placement methodology based on an objective function that minimizes the response reconstruction uncertainty at non-instrumented floors. The California Strong Motion Instrumentation Program (CSMIP) has been instrumenting civil engineering structures, including buildings and bridges, for over 50 years. One of the primary goals of this program is to provide high-quality structural response data to the engineering community. However, due to the high cost of sensor installation and maintenance, it is not feasible to densely instrument most structures. For example, many buildings are currently instrumented using engineering judgment, with only a few sensors distributed along their height. To make this process more quantitative and systematic, we developed a method for identifying optimal sensor locations along the height of building structures, enabling the interpolation of responses at non-instrumented floors with minimal uncertainty. The proposed method is limited to stationary excitations and buildings with uniform mass and stiffness distributions along their height. Under these assumptions, the kernel function of the Gaussian Process Regression (GPR) model can be analytically derived based on the shear-flexural beam model and random vibration theory. We verified and validated the proposed method using both simulated and real data. In addition, we proposed a table that can be used to support more quantitative and informed decision-making regarding sensor placement.