- Email: [email protected]
Multi-stage regression linear parametric models of room temperature in office buildings
Building and Environment 56 (2012) 69e77
Contents lists available at SciVerse ScienceDirect
Building and Environment journal homepage: www.elsevier.com/locate/buildenv
Multi-stage regression linear parametric models of room temperature in office buildings Siyu Wu, Jian-Qiao Sun* School of Engineering, University of California, 5200 North Lake Road, Merced, CA 95343, USA
a r t i c l e i n f o
a b s t r a c t
Article history: Received 25 October 2011 Received in revised form 3 January 2012 Accepted 22 February 2012
Mathematical models of the heating, ventilation, and air conditioning (HVAC) components play an important role in control design and fault detection of the system. The work in this paper incorporates architectural parameters in linear parametric models of room temperature in office buildings. Specifically, we allow the physics-based autoregression moving average (pbARMAX) model to have a multistage structure in order to explicitly include the architectural parameters of the room in the model. Extensive measurements of the room temperature are used to develop and validate the multi-stage model. The resulting model can predict the temperature in different rooms accurately in both shortterm and long-term. Over a period of four weeks, the predictions have a root mean squared error less than 0.10 with a coefficient of determination larger than 0.99. Ó 2012 Elsevier Ltd. All rights reserved.
Keywords: HVAC system Temperature model Multi-stage regression Autoregressive model
1. Introduction The significant energy consumed by heating, ventilation, and air conditioning (HVAC) systems [1e3] leads to a worldwide increase in relevant research over the past decades. According to the National Institute of Standards and Technology, various system faults, inefficient control strategies, and poor commissioning contribute 10%e40% share of the energy consumed by HVAC systems [4]. Mathematical models of HVAC units provide bases for detecting faults, upgrading control strategies, and improving commissioning, and thus have a potential to reduce energy consumption of HVAC systems by 20%e30% [5]. An important part of modeling HVAC systems is to develop accurate, robust, and yet simple models of the room temperature. The aim of this study is to develop such models to predict the room temperature in office buildings. Specifically, we present a modeling approach that allows the pbARMAX model to have a multi-stage structure in order to explicitly include the architectural parameters of the room in the model. The indoor temperature characterizes the thermal condition of a room, as is the case in the predicted mean vote (PMV) model [6,7], the ASHRAE Standard 55 [8] and the ISO Standard 7730 [9]. Ascione et al. need precise room temperature estimation to evaluate the strict thermo-hygrometric environment in museums [10]. Balaras
* Corresponding author. Tel.: þ1 209 228 4540. E-mail addresses: [email protected] (S. Wu), [email protected] (J.Q. Sun). 0360-1323/$ e see front matter Ó 2012 Elsevier Ltd. All rights reserved. doi:10.1016/j.buildenv.2012.02.026
et al. require the accurate knowledge of the room temperature to examine the indoor environment in hospitals [11]. The room temperature is also required for the analysis of thermal performance, indoor air quality [12], and for the evaluation of the indoor environment in naturally ventilated buildings [13]. Furthermore, a robust model that accurately predicts the room temperature is important for control design. Engdahl and Johansson have studied optimal supply air temperature in a variable air volume (VAV) system [14]. Yang and Kim predict the time of room temperature variations [15]. Orosa has investigated the thermal comfort based control strategy [16]. Tanimoto and Hagishima have developed a Markov model for the on-off cooling schedule control [17]. There are two types of room temperature models, i.e., analytical models and numerical models. The analytical models widely implemented in popular simulation tools such as EnergyPlus [18], DOE-2, TRANSYS, Modelica [19], etc. incorporate architectural parameters, and provide the most comprehensive description of the thermal process in the building with accurate estimation of various system outputs. The analytical models allow parametric studies of influence of the building envelope on indoor thermal behavior. This is very important because architectural and material parameters of a building greatly influence its thermal performance. Extensive parametric studies will help us optimize the design of and the material selection for buildings. However, various simplifying assumptions to deal with complexity of thermal interactions, unmeasured disturbances, uncertainty in thermal properties of structural elements and other parameters make it quite a challenge to obtain reliable analytical models [20].
S. Wu, J.Q. Sun / Building and Environment 56 (2012) 69e77
Nomenclature
DEh,d r q z ɸ
j e h k m _ m n q t x A, B, C AHU Cv Cp D H L
internal heat generated by human and devices (W) air density (kg m3) solar flux angle on windows ( ) window-to-wall ratio (%) heat gain from solar flux (W) solar heat flux density (W m2) prediction error heat convection coefficient (Wm2 K1) thermal conductivity coefficient (Wm1 K1) number of measurements mass flow rate (kg s1) time sequence of measurements delay operator time distance from inside surface (m) polynomials in the transfer functions for the ARMAX model air handling unit volumetric heat capacity of air (kJ kg1 K1) heat capacity of air at constant pressure (kJ kg1 K1) damper position (%) specific enthalpy of humid air (kJ kg1) thickness of wall or window (m)
On the contrary, data-driven, fully numerical modeling approaches establish models by using only input and output measurements. Many researchers have developed numerical models of the room temperature under various conditions. RíosMoreno et al. have demonstrated that the linear autoregression model with external input (ARX) model can be adopted to predict classroom indoor air temperature with very high accuracy [21]. Lowry and Lee have investigated the response of internal temperature in an office building and discovered that the output-error (OE) model provides reliable predictions [20]. Peitsman and Bakker have developed a multiple input single output ARX model to evaluate the performance of a VAV unit [22]. Yiu and Wang have studied the optimal order of a multiple input multiple output autogregressive moving average (ARMAX) model to forecast the performance of an air conditioning system of an office building in Hong Kong [23]. Mustafaraj et al. have compared different linear parametric models of the room temperature in an office, and find that the Box-Jenkins (BJ) model outperforms the ARX model and the ARMAX model under certain environmental conditions [24]. Although these numerical models are computationally efficient due to their simple structures, they are heavily dependent on the measurements, which implies poor generalization of the model in some parameter space. To help optimize the model selection, and to improve the generalization ability of fully numerical models, researchers have turned to semi-physical or gray-box modeling approach. The graybox models are constructed based on both the underlying physical laws and experimental data. Déqué and colleagues have implemented the gray-box technique to simulate the temperature variation in a ground floor flat [25]. Ghiaus et al. have conducted the gray-box identification of some temperature properties within the constant air AHU units [26]. Based on the physical laws, Wen and Smith have applied the gray-box approach to model the room temperature in VAVs [27]. Leephakpreeda has combined the gray prediction model and the adaptive comfort theory model to estimate the indoor comfort temperature [28]. Nevertheless, the gray-
P R S T V X RH VAV
parameter group residual of temperature simulation (K) surface area exposed to sun (m2) temperature (K) volume of room (m3) an attribute relative humidity of air (%) variable air volume unit
Subscripts i inside o outside oa outside air rm room wa wall wd window clg cooling water htg heating water surf surface disch discharge Superscript measurements e estimated values ˇ
70
box models are numerical and therefore don’t include architectural parameters. It is thus natural to consider combining the strength of the analytical model and the data-driven numerical approach. Tsilingiris has studied the influence of structural parameters of a wall on the heat loss [29]. Devgan et al. [30] and Aste et al. [31] have investigated the influence of external wall and window construction parameters on the building overall heat transfer and energy performance of well insulated buildings. By exploring the optimal dimensions of walls and glasses, Hwang and Shu [32], and Sozer [33] have tried to improve the thermal comfort, and the energy efficiency of buildings. Korolija et al. has concluded that it is not possible to form a reliable judgment about the building energy performance without considering the architectural parameters [34]. Numerical models without including the building parameters do not reveal the relationship between the building thermal performance and the parameters. We propose two new multi-stage regression, physics-based linear parametric (mpbARMAX) models of the room temperature to take the advantages of both analytical and numerical modeling approaches [35,36]. In other words, we build a two-stage regression structure into the mpbARMAX model to explicitly include building geometries. The rest of the paper is organized as follows. In Section 2, we review the pbARMAX model and introduce the two-stage regression structure. In Section 3, we develop the two-stage mpbARMAX models and validate them with extensive data from a building on the campus of University of California at Merced. The models are compared with a representative ARMAX model in both short-term and long-term predictions. Finally, Section 4 discusses and summarizes the findings. 2. Methodology We first review the two pbARMAX models of the room temperature. We then introduce linear regression representations of the coefficients in the pbARMAX models as a function of
S. Wu, J.Q. Sun / Building and Environment 56 (2012) 69e77
architectural parameters. The resulting mpbARMAX models have coefficients that are explicit functions of architectural parameters. 2.1. Physical equations The rate of energy change in a room can be expressed as [26,37,38],
dT _ p ðTdisch Trm Þ þ hwai Swa ðTwai Trm Þ ri VCv rm ¼ mC dt þ hwdi Swd Twdi Trm þ f þ DEh;d ;
Twao Twai xwa þ Twai ; Lwa
(2)
where Twao is the outside wall surface temperature, Lwa is the thickness of wall, and xwa is the distance from the inside wall surface. We further assume that in the thermal equilibrium, the heat convection on both the inside and outside surfaces of the wall is equal to the heat conduction through the wall [39]. This gives
dTwa kwa ¼ hwai ðTrm Twai Þ ¼ hwao ðTwao Toa Þ; dxwa
(3)
where kwa is the thermal conductivity coefficient of wall, and Toa is the outside air temperature. Eqs. (2) and (3) suggest that Twai can be expressed as a linear function of Trm and Toa, i.e.,
Twai ¼ P1 Trm þ P2 Toa ;
Twdi ¼ P3 Trm þ P4 Toa ;
(5)
where P3 and P4 are defined in Appendix. With these simplifications, Eq. (1) now reads,
dTrm ¼ Prm Trm þ Pdisch Tdisch þ Poa Toa þ DE; dt
(6)
where Prm, Pdisch, Poa and DE are defined in Appendix. 2.2. Numerical models of room temperature
(1)
where Trm, Tdisch, Twai and Twdi represent the temperature of the room, discharge air, inside wall surface, and inside window surface; V, Swa and Swd denote the room volume, surface area of the wall, and surface area of the window; ri, Cv, Cp, hwai and hwdi are discharge the air density, volumetric heat capacity of the air, and heat capacity of the air at constant pressure; f, and DEh,d represent the heat gain from solar flux, and from internal human and device. All the temperatures should be viewed as an average of the media. We don’t include the latent heat in this equation because there is a very low level of variations in humidity during the period of experiments reported later. However, the proposed modeling method can readily take the latent heat into consideration when _ p ðTdisch Trm Þ of the humidity data are available, i.e., replacing mC _ disch Hrm Þ, where H is the specific enthalpy of Eq. (1) with mðH humid air. Another simplification in the model is that the heat transfer through walls between rooms controlled by the same VAV with the same thermostat setpoint is neglected. This is based on the fact that, according to the data collected, the temperature difference between the room and the outside air is 2200% larger than the difference between two rooms on average. Therefore, the heat convection between the room and the outside air dominates the thermal interaction between the room and its surrounding. Nevertheless, our approach allows incorporating the heat convection between rooms into Eq. (1). In the following, we make several simplifying assumptions about the heat transfer through walls and windows. Take the wall as an example. We assume that the thermal conductivity of the wall is constant, and that the heat conduction through the wall is one dimensional in steady-state. Thus, the gradient of Twa along the thickness direction is constant. And the temperature in the wall Twa is therefore a linear function of the thickness coordinate,
Twa ¼
71
(4)
where the coefficients P1 and P2 are defined in Appendix. Following the same arguments, we obtain an expression for Twdi as
Let Dt be a sample interval and Trm (n) denote the sampled value of Trm (t) at time t ¼ nDt where n is a positive integer. Applying the central difference scheme to Eq. (6) leads to the pbARMAX model as
Trm ðn þ 1Þ ¼ Trm ðn 1Þ þ 2Dt½Prm Trm ðnÞ þ Pdisch Tdisch ðnÞ þ Poa Toa ðnÞ þ DEðnÞ:
(7)
This model has four independent coefficients to determine, i.e., Prm, Pdisch, Poa and DE. The model does not explicitly include all local _ This makes it difficult control variables such as the air flow rate m. to use the model for control design. We can improve the pbARMAX _ explicitly appears in model by rearranging the coefficients so that m the model. Referring to the definition of Prm in Appendix, we rewrite it as
_ p Prm ¼ P5 Swa hwai P2 þ hwdi zP4 P5 mC ¼ Poa Pdisch ;
(8)
and substitute Eq. (8) into Eq. (7). We have
Trm ðn þ 1Þ ¼ Trm ðn 1Þ þ 2Dt a$ðToa ðnÞ Trm ðnÞÞ
_ þ b$mðnÞðT disch ðnÞ Trm ðnÞÞ þ DEðnÞ ;
(9)
where a and b are in general functions of architectural parameters. Another way to arrange parameters can lead to a different form of the same model as
Trm ðn þ 1Þ ¼ Trm ðn 1Þ þ 2Dt½c$ðTdisch ðnÞ Trm ðnÞÞ þd$Toa ðnÞ þ DEðnÞ;
ð10Þ
where c and d are also functions of architectural parameters. For brevity, we shall refer to the models in Eqs. (9) and (10) as pbARMAX1 and pbARMAX2, respectively. We shall use a large set of measurements from VAVs involving a wide range of architectural parameters to determine the pbARMAX1,2 models. This results in a collection of coefficients a, b, c and d. 2.3. Regression model of the coefficients With the help of the expressions in Appendix, we can derive the following linear relationships for the coefficients a, b, c and d,
a ¼ p11 =V þ p12 $Swa =V þ p13 $Swd =V;
(11a)
b ¼ p21 =V þ p22 $Swa =V þ p23 $Swd =V;
(11b)
c ¼ p31 =V þ p32 $Swa =V þ p33 $Swd =V;
(11c)
d ¼ p41 =V þ p42 $Swa =V þ p43 $Swd =V;
(11d)
where V is the volume of the room, Swa is the area of the wall, Swd is the area of the windows, and pij (1 i 4, 1 j 3) are constants. It should be noted that these coefficients are also dependent on other geometrical and material parameters. In the present work, we shall focus on 1/V, Swa/V and Swd/V. In fact, the surface-area-to-volume
72
S. Wu, J.Q. Sun / Building and Environment 56 (2012) 69e77
ratios Swa/V and Swd/V are the driving forces to the thermodynamic processes that minimize the free energy [40,41]. Substituting these linear relationships into Eqs. (9) and (10), we obtain two-stage regression models of the room temperature, denoted as mpbARMAX1,2.
Trm ðnþ1Þ ¼ Trm ðn1Þþ2Dt ðp11 =V þp12 $Swa =V þp13 $Swd =VÞ$ðToa ðnÞTrm ðnÞÞþðp21 =V þp22 $Swa =V _ þp23 $Swd =VÞ$mðnÞðT disch ðnÞTrm ðnÞÞþ DEðnÞ ; (12) and
Trm ðn þ 1Þ ¼ Trm ðn 1Þ þ 2Dt½ðp31 =V þ p32 $Swa =V þ p33 $Swd =VÞ$ðTdisch ðnÞ Trm ðnÞÞ þ ðp41 =V þ p42 $Swa =V þ p43 $Swd =VÞ$Toa ðnÞ þ DEðnÞ:
(13)
It should be noted that the mpbARMAX1,2 models can be identified with one-stage regression by directly using Eqs. (12) and (13). However, the two-stage regression has its merits. (1) According to Baguley [42], Ahn [43] and Kristjansson et al. [44], it is beneficial to separate regression modeling of fixed and random effects. Here, the architectural parameters are fixed, and the temperature variations are random. The variability of temperature measurements from different VAVs is usually greater than the variability of the measurements from the same VAV. The regressions in two stages tackle these two variabilities separately and effectively [45]. (2) The number of independent data sets of 1/V, Swa/V and Swd/V is significantly smaller than that of the temperatures Tdisch (n), Toa (n) and Trm (n). By separating these variables in two stages of regression, the samples in each stage have similar distributions, thus improving the so-called compound symmetry of samples as discussed in [46,47]. This usually leads to better predictions. (3) Compared to the one-stage regression, the ratio of training samples to predictors in the two-stage regression increases more than 200%. This reduces the potential of overfitting [48,49]. In the following, we shall demonstrate how to determine the mpbARMAX1,2 models and their ability to predict the room temperature. 3. Case studies In this section, we demonstrate how to identify and validate the proposed mpbARMAX models. The comparison between the mpbARMAX models and a representative ARMAX model in the same domain will also be presented. 3.1. Building description Extensive past and current measurements of the HVAC system in the Science and Engineering (SE) building of UC Merced are available to us. The SE building is instrumented with a network of sensors and controls for the HVAC system. Such a highly instrumented building serves as an ideal living laboratory to support the research on energy efficiency. The SE building is a four-floor, southwest orientation, 19,666 gross square meter building with primary use as office and laboratory. Two heating and cooling water bridges, nine unit heaters and ten AHUs work together with sixtyone VAVs to control and regulate 374 rooms and spaces in the building. There are two major types of rooms, i.e., the faculty office and the conference room. The geometrical dimensions of
Fig. 1. The distribution of VAVs on the third floor of the SE building on UC merced campus.
the faculty office, and the conference room are 3.2004 4.4714 3.0480 (length width height) cubic meters, and 9.6012 5.1816 3.0480 cubic meters respectively. Fig. 1 shows the third floor map of the SE building with detailed distribution of VAV units. Table 1 presents the geometries of rooms regulated by VAVs under AHU 10 (A10). It should be noted that the volumes are calculated based on the room heights measured from the floor to the suspended ceiling. We select this height because the volume in between is the space where the heat exchange occurs. Also, some rooms have no surface exposed to the outside. Merced, located in the San Joaquin Valley of California, is in the Mediterranean Steppe ecoregion. Its climate exhibits rich variations during a year. The summer is hot and dry from June to August, and the winter is cold and rainy from November to April. The measurements of the HVAC system have rich dynamics. As for the local control systems, the chilling plant has a lead-lagstandby setup of three chillers, and the heating plant consists of three, dual fuel boilers. The pump speeds of both the plants are modulated by PID control loops to maintain differential pressures at the discharges. The AHU units apply PID controls to regulate the supply fan variable-frequency drive (VFD), return fan VFD, and supply air temperature. These control loops maintain the duct pressure, return fan discharge pressure, and supply air temperature at their setpoints. The VAV units implement two separate control loops, i.e., the cooling loop and the heating loop, to keep the temperature at setpoints. Both of the loops apply PI controls with three operational modes: heating, cooling and deadband. Table 1 The geometries of the rooms or spaces controlled by VAVs of AHU 10. VAV
Swa (m2)
Swd (m2)
V (m3)
VAV
Swa (m2)
Swd (m2)
V (m3)
V351 V352 V353 V354 V355 V356 V357 V359 V360 V361 V362 V363 V364 V365
11.9 3.0 8.1 13.5 1.4 6.7 0 13.6 8.1 19.0 19.0 8.8 4.8 0.4
19.0 15.5 24.3 9.4 4.7 4.7 0 6.8 48.6 10.2 10.2 30.4 24.5 46.7
78.7 78.7 75.1 73.9 181.9 50.2 97.4 132.4 244.7 130.8 130.8 101.0 151.5 151.5
V251 V252 V253 V254 V255 V256 V257 V258 V259 V260 V261 V262 V263 V264
7.6 1.2 5.2 8.2 2.0 0 0 13.6 2.4 19.0 19.0 8.8 4.8 0.4
50.5 20.2 24.3 9.4 1.8 0 0 6.8 35.9 10.2 10.2 30.4 24.5 46.7
151.6 80.8 43.6 58.1 208.7 94.9 61.9 132.4 165.0 130.8 130.8 101.0 151.5 151.5
S. Wu, J.Q. Sun / Building and Environment 56 (2012) 69e77
3.2. Data preprocessing Measurements from fifty-eight VAVs controlled by two AHUs are available for the development of the mpbARMAX models. From May 31 to September 16, 2010, measurements of Trm, Tdisch, Toa and _ of each VAV over 109 days have been collected with a sample rate m of 5 min. This results in 31,392 samples with 232 attributes. VAVs with intermittent sensor measurements are excluded. The local building energy management system (BEMS) applies different thermostat setpoints during the day, i.e., 5 a.m.e1 a.m. of the next day, and the night, i.e., 1 a.m.e5 a.m. Most of the occupants work during the day. There is no significant difference in occupancy between weekdays and weekends. We divide the data by day and night. This is consistent with the idea of temporal partition proposed in [50]. The samples are divided into a training set and a validation set. The identified regression model is verified by the validation set of the VAVs not included in the training set. The training set includes all VAVs controlled by AHU 9 (A9) and VAVs controlled by A10 on the third floor. The validation set includes VAVs controlled by A10 on the second floor. This yields 44 VAVs for training, and 14 VAVs for validation with 31 and 12 different combinations of the parameters 1/V, Swa/V and Swd/V in each set. 3.3. Model identification e stage 1 We consider several commonly used metrics of goodness-of-fit to quantify the modeling and prediction accuracy of the regression models including the root mean squared error (RMSE), mean absolute error (MAE), coefficient of determination (r2) and maximum absolute error (MaxAE) [49]. These metrics are defined as
" RMSE ¼
m 2 1 X ^ X ~ X i i m i¼1
#1=2 ;
(14a)
(14b)
h
P^ ~ P ^ P ~ i2 XiÞ m X Xi iXi ð
; r2 ¼ h P^ 2 P ^ 2 P ~2 P ~ 2 m Xi ð XiÞ m Xi Xi ^ ~ MaxAE ¼ maxX i X i ;
(14c)
(14d)
i
where Xi is a particular attribute in time domain.
RMSE of prediction (K)
Table 2 The coefficients of the pbARMAX1 model and prediction errors over four weeks during the day for 14 VAVs in the training set. VAV
a
b
RMSE
MAE
r2
MaxAE
V351 V352 V353 V354 V355 V356 V357 V359 V360 V361 V362 V363 V364 V365
0.0040 4.6931E-4 9.0462E-5 3.8409E-5 7.8922E-7 1.8222E-5 2.5264E-4 6.3169E-5 1.3155E-5 1.5193E-5 3.6004E-5 5.2251E-5 1.8862E-4 9.1219E-5
2.6179E-5 3.6583E-5 1.7770E-5 1.0488E-5 1.7027E-7 1.1515E-5 6.9368E-5 3.3502E-5 6.6516E-6 9.5484E-6 1.1201E-5 1.1045E-5 3.1082E-5 8.8934E-6
0.1442 0.1238 0.1377 0.1014 0.0599 0.1529 0.1391 0.1040 0.0812 0.0725 0.0888 0.1593 0.1285 0.1725
0.1010 0.0898 0.0859 0.0678 0.0358 0.0881 0.1090 0.0735 0.0532 0.0398 0.0586 0.0865 0.0934 0.1346
0.9797 0.9837 0.9811 0.9910 0.9877 0.9552 0.9683 0.9793 0.9904 0.9830 0.9939 0.9077 0.9816 0.9747
1.7669 0.8897 1.4507 1.0675 0.5005 1.4381 1.1035 1.2149 0.6809 1.3269 0.8215 2.5716 1.2174 1.7140
An optimal length of time duration exists for a given training set. We determine the optimal number of training samples for each VAV by changing the training length in the increment of days to find the smallest prediction error over the longest predictive period. Fig. 2 shows an example of how the optimal training length during the day for a VAV (V201) is determined. With these optimal training lengths, we determine the regression models in Eqs. (9) and (10) with the method of least squares, resulting in a set of model coefficients a, b, c and d. These model coefficients are listed in Tables 2 and 3. Fig. 3 shows an example of the modeling of the room temperature over two days by the pbARMAX1 model for V351. The modeling error is negligible with RMSE and r2 being 0.0057 and 0.9998 during the day, and 0.0378 and 0.9983 during the night. We have found that both the pbARMAX1 and pbARMAX2 models have a similar modeling accuracy for all the VAVs. 3.4. Model identification e stage 2
m 1 X ^ ~ ; MAE ¼ X i X i m i¼1
0.02
73
pbARMAX (left axis)
0.015
The first stage of regression produces 31 sets of model coefficients a, b, c and d. The second stage of regression finds the linear expressions in Eq. (11) that relate the coefficients a, b, c and d to the architectural parameters 1/V, Swa/V and Swd/V. The second stage continues to use VAVs in the training set. The method of least squares is used to determine the regression coefficients pij which are listed in Table 4. Table 5 lists the goodness-of-fit of the second stage regression. The average RMSE of all four model coefficients is as small as 1.0291E-4. For the VAVs in the training set, we conduct four-week predictions by using the mpbARMAX models in Eqs. (12) and (13). The Table 3 The coefficients of the pbARMAX2 model and prediction errors over four weeks during the day for 14 VAVs in the training set.
1
pbARMAX2 (right axis) 0.015
0.01
0.01
0.005
0.005 0
10
20 30 40 Length of training data (days)
50
0 60
Fig. 2. The root mean square error (RMSE) of prediction vs. the length of the training data. For V201, and most of the other fifty-seven VAVs, the optimal training length for the proposed pbARMAX models during the summer is between forty to fifty days.
VAV
c
d
RMSE
MAE
r2
MaxAE
V351 V352 V353 V354 V355 V356 V357 V359 V360 V361 V362 V363 V364 V365
2.2778E-4 0.0023 0.0011 7.1563E-4 0.0038 0.0033 6.4945E-4 7.4673E-4 0.0020 0.0011 8.7483E-5 4.3001E-4 0.0019 4.0561E-5
1.3034E-4 0.0019 8.7864E-4 5.7748E-4 0.0033 0.0028 5.3243E-4 6.1395E-4 0.0016 8.9718E-4 8.4688E-5 3.4604E-4 0.0016 1.1577E-5
0.1432 0.1088 0.1355 0.1005 0.0874 0.1092 0.0758 0.0801 0.0781 0.0677 0.0749 0.1486 0.0886 0.1565
0.1004 0.0706 0.0866 0.0670 0.0643 0.0696 0.0478 0.0512 0.0524 0.0384 0.0443 0.0791 0.0547 0.1214
0.9785 0.9867 0.9822 0.9910 0.9781 0.9770 0.9843 0.9837 0.9909 0.9844 0.9945 0.9145 0.9865 0.9774
2.0819 1.0351 1.4643 1.0565 0.5594 1.4521 0.9054 1.1594 0.6482 1.3188 0.7711 2.3908 0.9269 1.5981
S. Wu, J.Q. Sun / Building and Environment 56 (2012) 69e77
Room temperature (K)
Room temperature (K)
74
Day 298 296 15
16
17
16 Time (days)
17
298 Night 297
measurements pbARMAX1
15
Fig. 3. The pbARMAX1 model has small modeling errors for one example VAV during both the day and night. The RMSE and r2 are 0.0057 and 0.9998 during the day, and 0.0378 and 0.9983 during the night.
coefficients pij and the architectural parameters in Table 1 are used in the models. Take the mpbARMAX1 model as an example, Table 6 shows the goodness-of-fit of such predictions for four VAVs in the dean’s suite on the third floor, which have distinct combinations of the architectural parameters 1/V, Swa/V and Swd/V. We can also compare the predictions of the two-stage and one-stage regression models by examining Table 6 and Table 2 (Tables 2 and 3 also include the goodness-of-fits of long-term predictions of pbARMAX models which will be explained in Section 3.5). It is clear that the mpbARMAX1 model is slightly less accurate than the one-stage model. The mpbARMAX2 model demonstrates the similar performance in accurate long-term predictions. Table 4 The second stage regression coefficients identified from VAVs in the training set. pi,j i i i i
¼ ¼ ¼ ¼
1 2 3 4
j¼1
j¼2
j¼3
0.0283 0.1200 1.7832 1.4488
0.0005 0.0003 0.0329 0.0269
0.0008 0.0008 0.0097 0.0080
Table 5 Goodness-of-fit of the second stage regression of the model coefficients.
RMSE MaxAE
^ a
^ b
^c
^ d
3.6497E-5 1.0639E-4
5.6382E-6 8.9216E-6
2.0870E-4 3.6774E-4
1.7075E-4 3.0310E-4
We can now validate the mpbARMAX models with the valida^ ^c and d ^ ^; b; tion set. We form the models by using the coefficients a from Eq. (11) with only pij and architectural parameters, and use them to conduct both short-term and long-term predictions. As an example, Fig. 4 shows the one-day prediction of V255 with the mpbARMAX2 model. The prediction not only matches well the measured temperature T~ rm when it fluctuates within the thermostat setpoints, but also tracks nicely when the VAV brings the room temperature back within the setpoints. Goodness-of-fits of the four-week long-term predictions are presented in Table 7. For the four VAVs in the dean’s suite on the second floor, the average RMSE and r2 of the predictions are 0.2572 and 0.9568, respectively. Compared to the single-stage regression models, the validation errors of the two-stage regression models demonstrate no significant difference. For all fourteen VAVs in the validation set, the average RMSE and r2 of the predictions of the one-stage regression models are 15% more and 3% less than those of the two-stage regression models. These results show that the proposed models can indeed capture the underlying thermodynamics of the system.
3.6. Comparison with an ARMAX model According to Schunn and Wallach [49], one cannot apply absolute standards in assessing the quality of a particular goodness-offit value. Each goodness-of-fit should be compared to those obtained by previous models in the same domain. To this end, we shall compare the mpbARMAX1,2 models with existing models of the room temperature. Many linear parametric models such as ARX, ARMAX, BJ, OE, etc. have been successfully implemented to model the room temperature in HVAC units. In a recent work [24], an ARMAX model has been claimed to be very effective to predict the room temperature with relatively small errors in the summer and fall, the same periods when our data were collected. Following the modeling approach in [24], we test different combinations of na, nb, nc and nk, in which na, nb and nc range from 1 to 10 and nk ranges from 1 to 6. And we find that the ARMAX of order [na ¼ 1, nb ¼ 1, nc ¼ 1, nk ¼ 3], or ARMAX1113 for short, produces the best validation results on the SE building. This model is the same as one of the examples studied in [24]. The structure of the model ARMAX1113 is
Trm ðnÞ ¼
3.5. Model validation We first validate the pbARMAX1,2 models, which merely serves to assure that the appropriate models are chosen for the room temperature so that the mpbARMAX1,2 models are based on a solid foundation. We use the pbARMAX1,2 models to conduct four-week longterm predictions with the data from the training set. Take fourteen VAVs in the short wing on the third floor as examples. Tables 2and 3 in the previous Section 3.3 present the goodness-of-fits for the pbARMAX1,2 models. The average RMSE, MAE, r2, and MaxAE of the predictions are 0.12, 0.07, 0.97 and 1.26 for the pbARMAX1 model, while 0.10, 0.06, 0.98 and 1.24 for the pbARMAX2 model. Compared to the modeling errors presented in Section 3.3, the prediction demonstrates an increase in RMSE, and a slight decrease in r2. Nevertheless, the overall predictions of the pbARMAX1,2 models are accurate.
B1 ðqÞ B ðqÞ B ðqÞ _ nk Þ þ 2 Tdisch ðn nk Þ þ 3 Tclg ðn nk Þ mðn AðqÞ AðqÞ AðqÞ B4 ðqÞ B5 ðqÞ B ðqÞ þ Toa ðn nk Þ þ T ðn nk Þ þ 6 RHdisch AðqÞ AðqÞ htg AðqÞ B7 ðqÞ CðqÞ ðn nk Þ þ RHoa ðn nk Þ þ eðnÞ; ð15Þ AðqÞ AðqÞ
where Tclg and Thtg present the temperatures of cooling water, and heating water; RHdisch and RHoa denote the relative humidities of the discharge air and outside air; A(q), Bi(q) and C(q) are polynomials in q defined in Appendix. The variables in Eq. (15) that don’t appear in the mpbARMAX models are dropped in the comparison. Table 6 The four-week prediction errors generated by the mpbARMAX1 model for the example VAVs in the training set.
V351 V352 V353 V354
RMSE
MAE
r2
MaxAE
0.3541 0.1273 0.2591 0.4786
0.3204 0.0890 0.2330 0.4695
0.9752 0.9894 0.9888 0.9916
2.0186 1.2364 1.5983 1.3024
Room temperature (K)
S. Wu, J.Q. Sun / Building and Environment 56 (2012) 69e77
298
measurements training prediction
Training RMSE = 0.1098 r2 = 0.9303
297 296 Prediction RMSE = 0.1842 r2 = 0.9191
295 294 41
42
43
Time (days) Fig. 4. The short-term prediction of the two-stage regression mpbARMAX2 model over the validation set of VAVs demonstrates goodness-of-fits nearly as good as those over the training set.
Table 7 The four-week prediction errors generated by the mpbARMAX2 model for the example VAVs in the validation set.
V251 V252 V253 V254
RMSE
MAE
r2
MaxAE
0.3431 0.1462 0.2557 0.2839
0.3009 0.1035 0.1986 0.2755
0.9281 0.9800 0.9264 0.9930
1.4529 1.1708 2.1889 0.9272
The VAVs in the validation set are used for the comparison. The mpbARMAX1,2 models are formed as discussed previously. The ARMAX1113 model uses the historical data to train the model and obtains the optimal parameters. We apply the optimal training lengths for both approaches. We consider the prediction of the room temperature over a period of four weeks. Take V259 as an example. Fig. 5 compares the predictions from both models. The ARMAX1113 model performs well by providing a small RMSE as 0.0992, and a large r2 as 0.9712. The mpbARMAX1 model demonstrates an even better prediction with a 85.5% decrease in RMSE to 0.0143, and a 1.7% increase in r2 to 0.9878. Furthermore, the mpbARMAX1 model behaves better especially when the temperature is oscillating around the cooling setpoint. The observed MAE of the mpbARMAX1 model is 79.9% smaller than that of the ARMAX1113 model. For other thirteen VAVs in the validation set as well as the mpbARMAX2 model, we have observed the similar results. Table 8 compares the average goodness-of-fit of both the mpbARMAX1,2 model and the ARMAX1113 model for all VAVs in the
Fig. 5. The mpbARMAX1 model in the four-week long-term predictions over an example VAV outperforms the ARMAX1113 model in terms of all four validation metrics.
75
Table 8 The long-term prediction errors over four weeks by the mpbARMAX models (MR) and the ARMAX1113 model (AR). The mpbARMAX models demonstrate better performances consistently. MAE
r2
MaxAE
Prediction (weeks)
RMSE MR/AR
MR/AR
MR/AR
MR/AR
1 2 3 4
0.0420/0.0672 0.0541/0.1373 0.0591/0.1498 0.0908/0.1643 116%[/144%[
0.0382/0.0351 0.0384/0.0612 0.0431/0.0700 0.0806/0.0895 110%[/154%[
0.9995/0.9876 0.9980/0.9807 0.9974/0.9724 0.9959/0.9587 0.36%Y/3.92%Y
0.2414/0.6020 0.3335/1.1166 0.3984/1.1812 0.4853/1.7723 101%[/194%[
validation set. Compared to the ARMAX1113 model, the mpbARMAX1,2 models have a 19.5% smaller increase in RMSE over the four weeks, a 29.6% smaller increase in MAE, a 48.0% smaller increase in MaxAE, while a 90.9% less decrease in r2. This illustrates that the mpbARMAX1,2 models produce more accurate long-term predictions than the ARMAX1113 model. 4. Discussions and conclusions (1) Different training lengths in time for the first stage regression model lead to slightly different sets of coefficients a, b, c and d. The overall effect of such differences on the predictions of the first stage regression models is limited. Compared to the models identified with the optimal training lengths, the models identified with twenty days less or more than the optimal training length have a maximum 5% difference in goodness-of-fit of the prediction. (2) Our data were sampled at 5-min interval. All the results reported in this paper are based on the 5-min interval. Sampling time intervals, however, affect the outcome of the regression results. We have tested 1 min, 10 min, 15 min and 30 min in developing the mpbARMAX models. The shorter or longer sampling interval is achieved by linearly interpolating the original data or by deleting the extra data points. For sampling intervals longer than 5 min, the resulting mpbARMAX models have less than 5% change in the metrics for validation. The mpbARMAX models developed with the 5-min sampling interval on average produce the best predictions for all the VAVs. On the other hand, the mpbARMAX models developed with 1-min sampling interval produce the worse results compared to all other sampling intervals. (3) The proposed model is data-driven and can be updated when significant changes occur in HVAC control strategies, major equipment and operation profiles. One can even periodically update the model regardless, leading to an adaptive modeling strategy of HVAC systems. (4) The proposed model could also be used to predict other quantities of interest in HVAC systems such as the relative humidity. We limit the discussions to the temperature model in this paper only because the relative humidity measurements are unavailable. We are currently trying to manage a way to estimate the relative humidity based on the available measurements at the AHU level. We shall report this study in the future publications. To conclude, we have developed a new strategy to incorporate the architectural parameters into linear parametric regression models of the room temperature. The regression model is developed in two stages. In the first stage, linear parametric regression models of the room temperature based on the thermodynamic equations are developed with the data from VAVs involving different architectural parameters. The collection of the parameters
76
S. Wu, J.Q. Sun / Building and Environment 56 (2012) 69e77
of the regression models are then related to the architectural parameters with the help of another set of linear regression relationships. This is called the second stage regression. The model developed this way has been validated with the data in the validation set. We have found that the two-stage regression models are able to predict the room temperature in both short and long terms with a high accuracy. Furthermore, the proposed models demonstrate improved goodness-of-fit compared to a representative ARMAX model in the literature which has been proved to be very effective in the same application. Acknowledgment The authors would like to thank John Elliott of UC Merced for helping us to get the research started with access to the HVAC web control system. We deeply appreciate the generous technical guidance and input from John Elliott, Philip Haves and Michael Wetter of LBNL, and David Auslander of UC Berkeley. This paper was based on work of a CITRIS Seed Grant project. Appendix The parameters involved in the pbARMAX model are defined as follows.
! hwai Lwa k2wa 1þ ; P1 ¼ kwa þhwai Lwa kwa Lwa ðhwai þhwao ÞþL2wa hwai hwao P2 ¼
kwa hwao ; kwa ðhwai þ hwao Þ þ Lwa hwai hwao
hwdi Lwd P3 ¼ kwd þ hwdi Lwd
1þ
(16)
(17) !
k2wd
; kwd Lwd hwdi þ hwdo þ L2wd hwdi hwdo (18)
P4 ¼
P5 ¼
kwd hwdo ; kwd hwdi þ hwdo þ Lwd hwdi hwdo 1
ri VCv
(19)
;
(20)
_ p ; Prm ¼ P5 Swa hwai ðP1 1Þ þ hwdi zðP3 1Þ mC
(21)
_ p; Pdisch ¼ P5 mC
(22)
Poa ¼ P5 Swa hwai P2 þ hwdi zP4 ;
(23)
DE ¼ P5 f þ DEh;d :
(24)
The polynomials of the transfer functions in the ARMAX1113 model are given by
AðqÞ ¼ 1 þ a1 q1 þ a2 q2 þ . þ ana qna ; Bi ðqÞ ¼ 1 þ b1;i q1 þ b2;i q2 þ . þ bnb ;i qnb ;
(25) i ¼ 1; 2; .; 7 (26)
CðqÞ ¼ 1 þ c1 q1 þ c2 q2 þ . þ anc qnc ; where na, nb, and nc are the orders of the polynomials in q.
(27)
References [1] http://buildingsdatabook.eren.doe.gov/TableView.aspx?table¼1.1.1, U.S. residential and commercial buildings total primary energy consumption (quadrillion btu and percent of total). [2] http://buildingsdatabook.eren.doe.gov/TableView.aspx?table¼1.1.3, Buildings share of U.S. primary energy consumption (percent). [3] http://buildingsdatabook.eren.doe.gov/TableView.aspx?table¼1.1.9, Buildings share of U.S. electricity consumption (percent). [4] Schein J, Bushby ST, Castro NS, House JM. A rule-based fault detection method for air handling units. Energy Build 2006;38(12):1485e92. [5] Yoshida H, Kumar S. Rarx algorithm based model development and application to real time data for on-line fault detection in VAV AHU units. IBPSA Build Simulat 1999;99:161e8. [6] Fanger PO. Thermal comfort. Copenhagen: Danish Technical Press; 1970. [7] Fanger PO, Melikov AK, Hanzawa H, Ring J. Air turbulence and sensation of draught. Energy Build 1988;12(1):21e39. [8] American Society of Heating Refrigerating and Air Conditioning Engineers (ASHARE), Thermal environmental conditions for human occupancy (ASHARE, Standard 55-1992). [9] International Standards Organization (ISO), Moderate thermal environments: determination of the PMV and PPD indices and specification of the conditions for thermal comfort (ISO 7730:1994). [10] Ascione F, Bellia L, Capozzoli A, Minichiello F. Energy saving strategies in airconditioning for museums. Appl Thermal Eng 2009;29(4):676e86. [11] Balaras CA, Dascalaki E, Gaglia A. HVAC and indoor thermal conditions in hospital operating rooms. Energy Build 2007;39(4):454e70. [12] Yan D, Song F, Yang X, Jiang Y, Zhao B, Zhang X, et al. An integrated modeling tool for simultaneous analysis of thermal performance and indoor air quality in buildings. Build Environ 2008;43(3):287e93. [13] Wei S, Sun Y, Li M, Lin W, Zhao D, Shi Y, et al. Indoor thermal environment evaluations and parametric analyses in naturally ventilated buildings in dry season using a field survey and PMVe-PPDe model. Build Environ 2011;46(6): 1275e83. [14] Engdahl F, Johansson D. Optimal supply air temperature with respect to energy use in a variable air volume system. Energy Build 2004;36(3):205e18. [15] Yang IH, Kim KW. Prediction of the time of room air temperature descending for heating systems in buildings. Build Environ 2004;39(1):19e29. [16] Orosa JA. A new modelling methodology to control HVAC systems. Expert Syst Appl 2011;38(4):4505e13. [17] Tanimoto J, Hagishima A. State transition probability for the Markov model dealing with on/off cooling schedule in dwellings. Energy Build 2005;37(3): 181e7. [18] Pedrini A, Westphal FS, Lamberts R. A methodology for building energy modelling and calibration in warm climates. Build Environ 2002;37(8e9): 903e12. [19] Wetter M. Modelica-based modelling and simulation to support research and development in building energy and control systems. J Build Perform Simulation 2009;2(2):143e61. [20] Lowry G, Lee MW. Modelling the passive thermal response of a building using sparse BMS data. Appl Energy 2004;78(1):53e62. [21] Ríos-Moreno GJ, Trejo-Perea M, Castañeda-Miranda R, HernándezGuzmán VM, Herrera-Ruiz G. Modelling temperature in intelligent buildings by means of autoregressive models. Autom Construct 2007;16(5):713e22. [22] Peitsman HC, Bakker VE. Application of black-box models to HVAC systems for fault detection. ASHRAE Trans 1996;102(1):628e40. [23] Yiu JCM, Wang S. Multiple ARMAX modeling scheme for forecasting air conditioning system performance. Energy Convers Manag 2007;48(8): 2276e85. [24] Mustafaraj G, Chen J, Lowry G. Development of room temperature and relative humidity linear parametric models for an open office using BMS data. Energy Build 2010;42(3):348e56. [25] Déqué F, Ollivier F, Poblador A. Grey boxes used to represent buildings with a minimum number of geometric and thermal parameters. Energy Build 2000; 31(1):29e35. [26] Ghiaus C, Chicinas A, Inard C. Grey-box identification of air-handling unit elements. Control Eng Pract 2007;15:421e33. [27] Wen J, Smith TF. Development and validation of online models with parameter estimation for a building zone with VAV system. Energy Build 2007; 39(1):13e22. [28] Leephakpreeda T. Grey prediction on indoor comfort temperature for HVAC systems. Expert Syst Appl 2008;34(4):2284e9. [29] Tsilingiris PT. Wall heat loss from intermittently conditioned spaces e the dynamic influence of structural and operational parameters. Energy Build 2006;38(8):1022e31. [30] Devgan S, Jain AK, Bhattacharjee B. Predetermined overall thermal transfer value coefficients for composite, hot-dry and warm-humid climates. Energy Build 2010;42(10):1841e61. [31] Aste N, Angelotti A, Buzzetti M. The influence of the external walls thermal inertia on the energy performance of well insulated buildings. Energy Build 2009;41(11):1181e7. [32] Hwang RL, Shu SY. Building envelope regulations on thermal comfort in glass facade buildings and energy-saving potential for PMV-based comfort control. Build Environ 2011;46(4):824e34.
S. Wu, J.Q. Sun / Building and Environment 56 (2012) 69e77 [33] Sozer H. Improving energy efficiency through the design of the building envelope. Build Environ 2010;45(12):2581e93. [34] Korolija I, Marjanovic-Halburd L, Zhang Y, Hanby VI. Influence of building parameters and HVAC systems coupling on building energy performance. Energy Build 2011;43(6):1247e53. [35] Wang X, Eisenbrey J, Zeitz M, Sun JQ. Multi-stage regression analysis of acoustical properties of polyurethane foams. Sound Vib 2004;273(4e5): 1109e17. [36] Wang X, Sun JQ. On the fatigue analysis of non-Gaussian stress processes with asymmetric distribution. J Vib Acoust 2005;127(6):556e66. [37] Cunha JB, Couto C, Ruano AE. Real-time parameter estimation of dynamic temperature models for greenhouse environmental control. Control Eng Pract 1997;5(10):1473e81. [38] Jiménez MJ, Madsen H, Andersen KK. Identification of the main thermal characteristics of building components using MATLAB. Build Environ 2008; 43(2):170e80. [39] Bird RB, Stewart WE, Lightfoot EN. Transport phenomena. New York: Wiley; 2001. [40] Baierlein R. Thermal physics. Cambridge, England: Cambridge University Press; 2003. [41] Perrot P. A to Z of thermodynamics. Oxford, England: Oxford University Press; 1998.
77
[42] Baguley T. Beyond ANOVA: from repeated measures to multilevel models. In: Proceedings of PsyPAG mathematics statistics and computing workshop. Nottingham: Nottingham Trent University; 2008. [43] Ahn J. Beyond single equation regression analysis: path analysis and multistage regression analysis. Am J Pharm Educ 2002;66:37e42. [44] Kristjansson SD, Kircher JC, Webb AK. Multilevel models for repeated measures research designs in psychophysiology: an introduction to growth curve modeling. Psychophysiology 2007;44(5):728e36. [45] Bland JM, Altman DG. Statistics notes: correlations, regression, and repeated data. Br Med J 1994;308:896. [46] Mardia KV, Kent JT, Bibby JM. Multivariate analysis (probability and mathematical statistics). Waltham, Massachusetts: Academic Press; 1980. [47] Maxwell SE, Delaney HD. Designing experiments and analyzing data: a model comparison perspective. London: Routledge Academic; 2003. [48] Bishop CM. Neural networks for pattern recognition. Oxford: Clarendon Press; 1995. [49] Schunn CD, Wallach D. Evaluating goodness-of-fit in comparison of models to data. In: Tack W, editor. Psychologie der Kognition: Reden and vorträge anlässlich der emeritierung von Werner Tack. Saarbrueken, Germany: University of Saarland Press; 2005. p. 115e54. [50] Wu S, Sun JQ. Cross-level fault detection and diagnosis of building HVAC systems. Build Environ 2011;46(8):1558e66.
Contents lists available at SciVerse ScienceDirect
Building and Environment journal homepage: www.elsevier.com/locate/buildenv
Multi-stage regression linear parametric models of room temperature in office buildings Siyu Wu, Jian-Qiao Sun* School of Engineering, University of California, 5200 North Lake Road, Merced, CA 95343, USA
a r t i c l e i n f o
a b s t r a c t
Article history: Received 25 October 2011 Received in revised form 3 January 2012 Accepted 22 February 2012
Mathematical models of the heating, ventilation, and air conditioning (HVAC) components play an important role in control design and fault detection of the system. The work in this paper incorporates architectural parameters in linear parametric models of room temperature in office buildings. Specifically, we allow the physics-based autoregression moving average (pbARMAX) model to have a multistage structure in order to explicitly include the architectural parameters of the room in the model. Extensive measurements of the room temperature are used to develop and validate the multi-stage model. The resulting model can predict the temperature in different rooms accurately in both shortterm and long-term. Over a period of four weeks, the predictions have a root mean squared error less than 0.10 with a coefficient of determination larger than 0.99. Ó 2012 Elsevier Ltd. All rights reserved.
Keywords: HVAC system Temperature model Multi-stage regression Autoregressive model
1. Introduction The significant energy consumed by heating, ventilation, and air conditioning (HVAC) systems [1e3] leads to a worldwide increase in relevant research over the past decades. According to the National Institute of Standards and Technology, various system faults, inefficient control strategies, and poor commissioning contribute 10%e40% share of the energy consumed by HVAC systems [4]. Mathematical models of HVAC units provide bases for detecting faults, upgrading control strategies, and improving commissioning, and thus have a potential to reduce energy consumption of HVAC systems by 20%e30% [5]. An important part of modeling HVAC systems is to develop accurate, robust, and yet simple models of the room temperature. The aim of this study is to develop such models to predict the room temperature in office buildings. Specifically, we present a modeling approach that allows the pbARMAX model to have a multi-stage structure in order to explicitly include the architectural parameters of the room in the model. The indoor temperature characterizes the thermal condition of a room, as is the case in the predicted mean vote (PMV) model [6,7], the ASHRAE Standard 55 [8] and the ISO Standard 7730 [9]. Ascione et al. need precise room temperature estimation to evaluate the strict thermo-hygrometric environment in museums [10]. Balaras
* Corresponding author. Tel.: þ1 209 228 4540. E-mail addresses: [email protected] (S. Wu), [email protected] (J.Q. Sun). 0360-1323/$ e see front matter Ó 2012 Elsevier Ltd. All rights reserved. doi:10.1016/j.buildenv.2012.02.026
et al. require the accurate knowledge of the room temperature to examine the indoor environment in hospitals [11]. The room temperature is also required for the analysis of thermal performance, indoor air quality [12], and for the evaluation of the indoor environment in naturally ventilated buildings [13]. Furthermore, a robust model that accurately predicts the room temperature is important for control design. Engdahl and Johansson have studied optimal supply air temperature in a variable air volume (VAV) system [14]. Yang and Kim predict the time of room temperature variations [15]. Orosa has investigated the thermal comfort based control strategy [16]. Tanimoto and Hagishima have developed a Markov model for the on-off cooling schedule control [17]. There are two types of room temperature models, i.e., analytical models and numerical models. The analytical models widely implemented in popular simulation tools such as EnergyPlus [18], DOE-2, TRANSYS, Modelica [19], etc. incorporate architectural parameters, and provide the most comprehensive description of the thermal process in the building with accurate estimation of various system outputs. The analytical models allow parametric studies of influence of the building envelope on indoor thermal behavior. This is very important because architectural and material parameters of a building greatly influence its thermal performance. Extensive parametric studies will help us optimize the design of and the material selection for buildings. However, various simplifying assumptions to deal with complexity of thermal interactions, unmeasured disturbances, uncertainty in thermal properties of structural elements and other parameters make it quite a challenge to obtain reliable analytical models [20].
S. Wu, J.Q. Sun / Building and Environment 56 (2012) 69e77
Nomenclature
DEh,d r q z ɸ
j e h k m _ m n q t x A, B, C AHU Cv Cp D H L
internal heat generated by human and devices (W) air density (kg m3) solar flux angle on windows ( ) window-to-wall ratio (%) heat gain from solar flux (W) solar heat flux density (W m2) prediction error heat convection coefficient (Wm2 K1) thermal conductivity coefficient (Wm1 K1) number of measurements mass flow rate (kg s1) time sequence of measurements delay operator time distance from inside surface (m) polynomials in the transfer functions for the ARMAX model air handling unit volumetric heat capacity of air (kJ kg1 K1) heat capacity of air at constant pressure (kJ kg1 K1) damper position (%) specific enthalpy of humid air (kJ kg1) thickness of wall or window (m)
On the contrary, data-driven, fully numerical modeling approaches establish models by using only input and output measurements. Many researchers have developed numerical models of the room temperature under various conditions. RíosMoreno et al. have demonstrated that the linear autoregression model with external input (ARX) model can be adopted to predict classroom indoor air temperature with very high accuracy [21]. Lowry and Lee have investigated the response of internal temperature in an office building and discovered that the output-error (OE) model provides reliable predictions [20]. Peitsman and Bakker have developed a multiple input single output ARX model to evaluate the performance of a VAV unit [22]. Yiu and Wang have studied the optimal order of a multiple input multiple output autogregressive moving average (ARMAX) model to forecast the performance of an air conditioning system of an office building in Hong Kong [23]. Mustafaraj et al. have compared different linear parametric models of the room temperature in an office, and find that the Box-Jenkins (BJ) model outperforms the ARX model and the ARMAX model under certain environmental conditions [24]. Although these numerical models are computationally efficient due to their simple structures, they are heavily dependent on the measurements, which implies poor generalization of the model in some parameter space. To help optimize the model selection, and to improve the generalization ability of fully numerical models, researchers have turned to semi-physical or gray-box modeling approach. The graybox models are constructed based on both the underlying physical laws and experimental data. Déqué and colleagues have implemented the gray-box technique to simulate the temperature variation in a ground floor flat [25]. Ghiaus et al. have conducted the gray-box identification of some temperature properties within the constant air AHU units [26]. Based on the physical laws, Wen and Smith have applied the gray-box approach to model the room temperature in VAVs [27]. Leephakpreeda has combined the gray prediction model and the adaptive comfort theory model to estimate the indoor comfort temperature [28]. Nevertheless, the gray-
P R S T V X RH VAV
parameter group residual of temperature simulation (K) surface area exposed to sun (m2) temperature (K) volume of room (m3) an attribute relative humidity of air (%) variable air volume unit
Subscripts i inside o outside oa outside air rm room wa wall wd window clg cooling water htg heating water surf surface disch discharge Superscript measurements e estimated values ˇ
70
box models are numerical and therefore don’t include architectural parameters. It is thus natural to consider combining the strength of the analytical model and the data-driven numerical approach. Tsilingiris has studied the influence of structural parameters of a wall on the heat loss [29]. Devgan et al. [30] and Aste et al. [31] have investigated the influence of external wall and window construction parameters on the building overall heat transfer and energy performance of well insulated buildings. By exploring the optimal dimensions of walls and glasses, Hwang and Shu [32], and Sozer [33] have tried to improve the thermal comfort, and the energy efficiency of buildings. Korolija et al. has concluded that it is not possible to form a reliable judgment about the building energy performance without considering the architectural parameters [34]. Numerical models without including the building parameters do not reveal the relationship between the building thermal performance and the parameters. We propose two new multi-stage regression, physics-based linear parametric (mpbARMAX) models of the room temperature to take the advantages of both analytical and numerical modeling approaches [35,36]. In other words, we build a two-stage regression structure into the mpbARMAX model to explicitly include building geometries. The rest of the paper is organized as follows. In Section 2, we review the pbARMAX model and introduce the two-stage regression structure. In Section 3, we develop the two-stage mpbARMAX models and validate them with extensive data from a building on the campus of University of California at Merced. The models are compared with a representative ARMAX model in both short-term and long-term predictions. Finally, Section 4 discusses and summarizes the findings. 2. Methodology We first review the two pbARMAX models of the room temperature. We then introduce linear regression representations of the coefficients in the pbARMAX models as a function of
S. Wu, J.Q. Sun / Building and Environment 56 (2012) 69e77
architectural parameters. The resulting mpbARMAX models have coefficients that are explicit functions of architectural parameters. 2.1. Physical equations The rate of energy change in a room can be expressed as [26,37,38],
dT _ p ðTdisch Trm Þ þ hwai Swa ðTwai Trm Þ ri VCv rm ¼ mC dt þ hwdi Swd Twdi Trm þ f þ DEh;d ;
Twao Twai xwa þ Twai ; Lwa
(2)
where Twao is the outside wall surface temperature, Lwa is the thickness of wall, and xwa is the distance from the inside wall surface. We further assume that in the thermal equilibrium, the heat convection on both the inside and outside surfaces of the wall is equal to the heat conduction through the wall [39]. This gives
dTwa kwa ¼ hwai ðTrm Twai Þ ¼ hwao ðTwao Toa Þ; dxwa
(3)
where kwa is the thermal conductivity coefficient of wall, and Toa is the outside air temperature. Eqs. (2) and (3) suggest that Twai can be expressed as a linear function of Trm and Toa, i.e.,
Twai ¼ P1 Trm þ P2 Toa ;
Twdi ¼ P3 Trm þ P4 Toa ;
(5)
where P3 and P4 are defined in Appendix. With these simplifications, Eq. (1) now reads,
dTrm ¼ Prm Trm þ Pdisch Tdisch þ Poa Toa þ DE; dt
(6)
where Prm, Pdisch, Poa and DE are defined in Appendix. 2.2. Numerical models of room temperature
(1)
where Trm, Tdisch, Twai and Twdi represent the temperature of the room, discharge air, inside wall surface, and inside window surface; V, Swa and Swd denote the room volume, surface area of the wall, and surface area of the window; ri, Cv, Cp, hwai and hwdi are discharge the air density, volumetric heat capacity of the air, and heat capacity of the air at constant pressure; f, and DEh,d represent the heat gain from solar flux, and from internal human and device. All the temperatures should be viewed as an average of the media. We don’t include the latent heat in this equation because there is a very low level of variations in humidity during the period of experiments reported later. However, the proposed modeling method can readily take the latent heat into consideration when _ p ðTdisch Trm Þ of the humidity data are available, i.e., replacing mC _ disch Hrm Þ, where H is the specific enthalpy of Eq. (1) with mðH humid air. Another simplification in the model is that the heat transfer through walls between rooms controlled by the same VAV with the same thermostat setpoint is neglected. This is based on the fact that, according to the data collected, the temperature difference between the room and the outside air is 2200% larger than the difference between two rooms on average. Therefore, the heat convection between the room and the outside air dominates the thermal interaction between the room and its surrounding. Nevertheless, our approach allows incorporating the heat convection between rooms into Eq. (1). In the following, we make several simplifying assumptions about the heat transfer through walls and windows. Take the wall as an example. We assume that the thermal conductivity of the wall is constant, and that the heat conduction through the wall is one dimensional in steady-state. Thus, the gradient of Twa along the thickness direction is constant. And the temperature in the wall Twa is therefore a linear function of the thickness coordinate,
Twa ¼
71
(4)
where the coefficients P1 and P2 are defined in Appendix. Following the same arguments, we obtain an expression for Twdi as
Let Dt be a sample interval and Trm (n) denote the sampled value of Trm (t) at time t ¼ nDt where n is a positive integer. Applying the central difference scheme to Eq. (6) leads to the pbARMAX model as
Trm ðn þ 1Þ ¼ Trm ðn 1Þ þ 2Dt½Prm Trm ðnÞ þ Pdisch Tdisch ðnÞ þ Poa Toa ðnÞ þ DEðnÞ:
(7)
This model has four independent coefficients to determine, i.e., Prm, Pdisch, Poa and DE. The model does not explicitly include all local _ This makes it difficult control variables such as the air flow rate m. to use the model for control design. We can improve the pbARMAX _ explicitly appears in model by rearranging the coefficients so that m the model. Referring to the definition of Prm in Appendix, we rewrite it as
_ p Prm ¼ P5 Swa hwai P2 þ hwdi zP4 P5 mC ¼ Poa Pdisch ;
(8)
and substitute Eq. (8) into Eq. (7). We have
Trm ðn þ 1Þ ¼ Trm ðn 1Þ þ 2Dt a$ðToa ðnÞ Trm ðnÞÞ
_ þ b$mðnÞðT disch ðnÞ Trm ðnÞÞ þ DEðnÞ ;
(9)
where a and b are in general functions of architectural parameters. Another way to arrange parameters can lead to a different form of the same model as
Trm ðn þ 1Þ ¼ Trm ðn 1Þ þ 2Dt½c$ðTdisch ðnÞ Trm ðnÞÞ þd$Toa ðnÞ þ DEðnÞ;
ð10Þ
where c and d are also functions of architectural parameters. For brevity, we shall refer to the models in Eqs. (9) and (10) as pbARMAX1 and pbARMAX2, respectively. We shall use a large set of measurements from VAVs involving a wide range of architectural parameters to determine the pbARMAX1,2 models. This results in a collection of coefficients a, b, c and d. 2.3. Regression model of the coefficients With the help of the expressions in Appendix, we can derive the following linear relationships for the coefficients a, b, c and d,
a ¼ p11 =V þ p12 $Swa =V þ p13 $Swd =V;
(11a)
b ¼ p21 =V þ p22 $Swa =V þ p23 $Swd =V;
(11b)
c ¼ p31 =V þ p32 $Swa =V þ p33 $Swd =V;
(11c)
d ¼ p41 =V þ p42 $Swa =V þ p43 $Swd =V;
(11d)
where V is the volume of the room, Swa is the area of the wall, Swd is the area of the windows, and pij (1 i 4, 1 j 3) are constants. It should be noted that these coefficients are also dependent on other geometrical and material parameters. In the present work, we shall focus on 1/V, Swa/V and Swd/V. In fact, the surface-area-to-volume
72
S. Wu, J.Q. Sun / Building and Environment 56 (2012) 69e77
ratios Swa/V and Swd/V are the driving forces to the thermodynamic processes that minimize the free energy [40,41]. Substituting these linear relationships into Eqs. (9) and (10), we obtain two-stage regression models of the room temperature, denoted as mpbARMAX1,2.
Trm ðnþ1Þ ¼ Trm ðn1Þþ2Dt ðp11 =V þp12 $Swa =V þp13 $Swd =VÞ$ðToa ðnÞTrm ðnÞÞþðp21 =V þp22 $Swa =V _ þp23 $Swd =VÞ$mðnÞðT disch ðnÞTrm ðnÞÞþ DEðnÞ ; (12) and
Trm ðn þ 1Þ ¼ Trm ðn 1Þ þ 2Dt½ðp31 =V þ p32 $Swa =V þ p33 $Swd =VÞ$ðTdisch ðnÞ Trm ðnÞÞ þ ðp41 =V þ p42 $Swa =V þ p43 $Swd =VÞ$Toa ðnÞ þ DEðnÞ:
(13)
It should be noted that the mpbARMAX1,2 models can be identified with one-stage regression by directly using Eqs. (12) and (13). However, the two-stage regression has its merits. (1) According to Baguley [42], Ahn [43] and Kristjansson et al. [44], it is beneficial to separate regression modeling of fixed and random effects. Here, the architectural parameters are fixed, and the temperature variations are random. The variability of temperature measurements from different VAVs is usually greater than the variability of the measurements from the same VAV. The regressions in two stages tackle these two variabilities separately and effectively [45]. (2) The number of independent data sets of 1/V, Swa/V and Swd/V is significantly smaller than that of the temperatures Tdisch (n), Toa (n) and Trm (n). By separating these variables in two stages of regression, the samples in each stage have similar distributions, thus improving the so-called compound symmetry of samples as discussed in [46,47]. This usually leads to better predictions. (3) Compared to the one-stage regression, the ratio of training samples to predictors in the two-stage regression increases more than 200%. This reduces the potential of overfitting [48,49]. In the following, we shall demonstrate how to determine the mpbARMAX1,2 models and their ability to predict the room temperature. 3. Case studies In this section, we demonstrate how to identify and validate the proposed mpbARMAX models. The comparison between the mpbARMAX models and a representative ARMAX model in the same domain will also be presented. 3.1. Building description Extensive past and current measurements of the HVAC system in the Science and Engineering (SE) building of UC Merced are available to us. The SE building is instrumented with a network of sensors and controls for the HVAC system. Such a highly instrumented building serves as an ideal living laboratory to support the research on energy efficiency. The SE building is a four-floor, southwest orientation, 19,666 gross square meter building with primary use as office and laboratory. Two heating and cooling water bridges, nine unit heaters and ten AHUs work together with sixtyone VAVs to control and regulate 374 rooms and spaces in the building. There are two major types of rooms, i.e., the faculty office and the conference room. The geometrical dimensions of
Fig. 1. The distribution of VAVs on the third floor of the SE building on UC merced campus.
the faculty office, and the conference room are 3.2004 4.4714 3.0480 (length width height) cubic meters, and 9.6012 5.1816 3.0480 cubic meters respectively. Fig. 1 shows the third floor map of the SE building with detailed distribution of VAV units. Table 1 presents the geometries of rooms regulated by VAVs under AHU 10 (A10). It should be noted that the volumes are calculated based on the room heights measured from the floor to the suspended ceiling. We select this height because the volume in between is the space where the heat exchange occurs. Also, some rooms have no surface exposed to the outside. Merced, located in the San Joaquin Valley of California, is in the Mediterranean Steppe ecoregion. Its climate exhibits rich variations during a year. The summer is hot and dry from June to August, and the winter is cold and rainy from November to April. The measurements of the HVAC system have rich dynamics. As for the local control systems, the chilling plant has a lead-lagstandby setup of three chillers, and the heating plant consists of three, dual fuel boilers. The pump speeds of both the plants are modulated by PID control loops to maintain differential pressures at the discharges. The AHU units apply PID controls to regulate the supply fan variable-frequency drive (VFD), return fan VFD, and supply air temperature. These control loops maintain the duct pressure, return fan discharge pressure, and supply air temperature at their setpoints. The VAV units implement two separate control loops, i.e., the cooling loop and the heating loop, to keep the temperature at setpoints. Both of the loops apply PI controls with three operational modes: heating, cooling and deadband. Table 1 The geometries of the rooms or spaces controlled by VAVs of AHU 10. VAV
Swa (m2)
Swd (m2)
V (m3)
VAV
Swa (m2)
Swd (m2)
V (m3)
V351 V352 V353 V354 V355 V356 V357 V359 V360 V361 V362 V363 V364 V365
11.9 3.0 8.1 13.5 1.4 6.7 0 13.6 8.1 19.0 19.0 8.8 4.8 0.4
19.0 15.5 24.3 9.4 4.7 4.7 0 6.8 48.6 10.2 10.2 30.4 24.5 46.7
78.7 78.7 75.1 73.9 181.9 50.2 97.4 132.4 244.7 130.8 130.8 101.0 151.5 151.5
V251 V252 V253 V254 V255 V256 V257 V258 V259 V260 V261 V262 V263 V264
7.6 1.2 5.2 8.2 2.0 0 0 13.6 2.4 19.0 19.0 8.8 4.8 0.4
50.5 20.2 24.3 9.4 1.8 0 0 6.8 35.9 10.2 10.2 30.4 24.5 46.7
151.6 80.8 43.6 58.1 208.7 94.9 61.9 132.4 165.0 130.8 130.8 101.0 151.5 151.5
S. Wu, J.Q. Sun / Building and Environment 56 (2012) 69e77
3.2. Data preprocessing Measurements from fifty-eight VAVs controlled by two AHUs are available for the development of the mpbARMAX models. From May 31 to September 16, 2010, measurements of Trm, Tdisch, Toa and _ of each VAV over 109 days have been collected with a sample rate m of 5 min. This results in 31,392 samples with 232 attributes. VAVs with intermittent sensor measurements are excluded. The local building energy management system (BEMS) applies different thermostat setpoints during the day, i.e., 5 a.m.e1 a.m. of the next day, and the night, i.e., 1 a.m.e5 a.m. Most of the occupants work during the day. There is no significant difference in occupancy between weekdays and weekends. We divide the data by day and night. This is consistent with the idea of temporal partition proposed in [50]. The samples are divided into a training set and a validation set. The identified regression model is verified by the validation set of the VAVs not included in the training set. The training set includes all VAVs controlled by AHU 9 (A9) and VAVs controlled by A10 on the third floor. The validation set includes VAVs controlled by A10 on the second floor. This yields 44 VAVs for training, and 14 VAVs for validation with 31 and 12 different combinations of the parameters 1/V, Swa/V and Swd/V in each set. 3.3. Model identification e stage 1 We consider several commonly used metrics of goodness-of-fit to quantify the modeling and prediction accuracy of the regression models including the root mean squared error (RMSE), mean absolute error (MAE), coefficient of determination (r2) and maximum absolute error (MaxAE) [49]. These metrics are defined as
" RMSE ¼
m 2 1 X ^ X ~ X i i m i¼1
#1=2 ;
(14a)
(14b)
h
P^ ~ P ^ P ~ i2 XiÞ m X Xi iXi ð
; r2 ¼ h P^ 2 P ^ 2 P ~2 P ~ 2 m Xi ð XiÞ m Xi Xi ^ ~ MaxAE ¼ maxX i X i ;
(14c)
(14d)
i
where Xi is a particular attribute in time domain.
RMSE of prediction (K)
Table 2 The coefficients of the pbARMAX1 model and prediction errors over four weeks during the day for 14 VAVs in the training set. VAV
a
b
RMSE
MAE
r2
MaxAE
V351 V352 V353 V354 V355 V356 V357 V359 V360 V361 V362 V363 V364 V365
0.0040 4.6931E-4 9.0462E-5 3.8409E-5 7.8922E-7 1.8222E-5 2.5264E-4 6.3169E-5 1.3155E-5 1.5193E-5 3.6004E-5 5.2251E-5 1.8862E-4 9.1219E-5
2.6179E-5 3.6583E-5 1.7770E-5 1.0488E-5 1.7027E-7 1.1515E-5 6.9368E-5 3.3502E-5 6.6516E-6 9.5484E-6 1.1201E-5 1.1045E-5 3.1082E-5 8.8934E-6
0.1442 0.1238 0.1377 0.1014 0.0599 0.1529 0.1391 0.1040 0.0812 0.0725 0.0888 0.1593 0.1285 0.1725
0.1010 0.0898 0.0859 0.0678 0.0358 0.0881 0.1090 0.0735 0.0532 0.0398 0.0586 0.0865 0.0934 0.1346
0.9797 0.9837 0.9811 0.9910 0.9877 0.9552 0.9683 0.9793 0.9904 0.9830 0.9939 0.9077 0.9816 0.9747
1.7669 0.8897 1.4507 1.0675 0.5005 1.4381 1.1035 1.2149 0.6809 1.3269 0.8215 2.5716 1.2174 1.7140
An optimal length of time duration exists for a given training set. We determine the optimal number of training samples for each VAV by changing the training length in the increment of days to find the smallest prediction error over the longest predictive period. Fig. 2 shows an example of how the optimal training length during the day for a VAV (V201) is determined. With these optimal training lengths, we determine the regression models in Eqs. (9) and (10) with the method of least squares, resulting in a set of model coefficients a, b, c and d. These model coefficients are listed in Tables 2 and 3. Fig. 3 shows an example of the modeling of the room temperature over two days by the pbARMAX1 model for V351. The modeling error is negligible with RMSE and r2 being 0.0057 and 0.9998 during the day, and 0.0378 and 0.9983 during the night. We have found that both the pbARMAX1 and pbARMAX2 models have a similar modeling accuracy for all the VAVs. 3.4. Model identification e stage 2
m 1 X ^ ~ ; MAE ¼ X i X i m i¼1
0.02
73
pbARMAX (left axis)
0.015
The first stage of regression produces 31 sets of model coefficients a, b, c and d. The second stage of regression finds the linear expressions in Eq. (11) that relate the coefficients a, b, c and d to the architectural parameters 1/V, Swa/V and Swd/V. The second stage continues to use VAVs in the training set. The method of least squares is used to determine the regression coefficients pij which are listed in Table 4. Table 5 lists the goodness-of-fit of the second stage regression. The average RMSE of all four model coefficients is as small as 1.0291E-4. For the VAVs in the training set, we conduct four-week predictions by using the mpbARMAX models in Eqs. (12) and (13). The Table 3 The coefficients of the pbARMAX2 model and prediction errors over four weeks during the day for 14 VAVs in the training set.
1
pbARMAX2 (right axis) 0.015
0.01
0.01
0.005
0.005 0
10
20 30 40 Length of training data (days)
50
0 60
Fig. 2. The root mean square error (RMSE) of prediction vs. the length of the training data. For V201, and most of the other fifty-seven VAVs, the optimal training length for the proposed pbARMAX models during the summer is between forty to fifty days.
VAV
c
d
RMSE
MAE
r2
MaxAE
V351 V352 V353 V354 V355 V356 V357 V359 V360 V361 V362 V363 V364 V365
2.2778E-4 0.0023 0.0011 7.1563E-4 0.0038 0.0033 6.4945E-4 7.4673E-4 0.0020 0.0011 8.7483E-5 4.3001E-4 0.0019 4.0561E-5
1.3034E-4 0.0019 8.7864E-4 5.7748E-4 0.0033 0.0028 5.3243E-4 6.1395E-4 0.0016 8.9718E-4 8.4688E-5 3.4604E-4 0.0016 1.1577E-5
0.1432 0.1088 0.1355 0.1005 0.0874 0.1092 0.0758 0.0801 0.0781 0.0677 0.0749 0.1486 0.0886 0.1565
0.1004 0.0706 0.0866 0.0670 0.0643 0.0696 0.0478 0.0512 0.0524 0.0384 0.0443 0.0791 0.0547 0.1214
0.9785 0.9867 0.9822 0.9910 0.9781 0.9770 0.9843 0.9837 0.9909 0.9844 0.9945 0.9145 0.9865 0.9774
2.0819 1.0351 1.4643 1.0565 0.5594 1.4521 0.9054 1.1594 0.6482 1.3188 0.7711 2.3908 0.9269 1.5981
S. Wu, J.Q. Sun / Building and Environment 56 (2012) 69e77
Room temperature (K)
Room temperature (K)
74
Day 298 296 15
16
17
16 Time (days)
17
298 Night 297
measurements pbARMAX1
15
Fig. 3. The pbARMAX1 model has small modeling errors for one example VAV during both the day and night. The RMSE and r2 are 0.0057 and 0.9998 during the day, and 0.0378 and 0.9983 during the night.
coefficients pij and the architectural parameters in Table 1 are used in the models. Take the mpbARMAX1 model as an example, Table 6 shows the goodness-of-fit of such predictions for four VAVs in the dean’s suite on the third floor, which have distinct combinations of the architectural parameters 1/V, Swa/V and Swd/V. We can also compare the predictions of the two-stage and one-stage regression models by examining Table 6 and Table 2 (Tables 2 and 3 also include the goodness-of-fits of long-term predictions of pbARMAX models which will be explained in Section 3.5). It is clear that the mpbARMAX1 model is slightly less accurate than the one-stage model. The mpbARMAX2 model demonstrates the similar performance in accurate long-term predictions. Table 4 The second stage regression coefficients identified from VAVs in the training set. pi,j i i i i
¼ ¼ ¼ ¼
1 2 3 4
j¼1
j¼2
j¼3
0.0283 0.1200 1.7832 1.4488
0.0005 0.0003 0.0329 0.0269
0.0008 0.0008 0.0097 0.0080
Table 5 Goodness-of-fit of the second stage regression of the model coefficients.
RMSE MaxAE
^ a
^ b
^c
^ d
3.6497E-5 1.0639E-4
5.6382E-6 8.9216E-6
2.0870E-4 3.6774E-4
1.7075E-4 3.0310E-4
We can now validate the mpbARMAX models with the valida^ ^c and d ^ ^; b; tion set. We form the models by using the coefficients a from Eq. (11) with only pij and architectural parameters, and use them to conduct both short-term and long-term predictions. As an example, Fig. 4 shows the one-day prediction of V255 with the mpbARMAX2 model. The prediction not only matches well the measured temperature T~ rm when it fluctuates within the thermostat setpoints, but also tracks nicely when the VAV brings the room temperature back within the setpoints. Goodness-of-fits of the four-week long-term predictions are presented in Table 7. For the four VAVs in the dean’s suite on the second floor, the average RMSE and r2 of the predictions are 0.2572 and 0.9568, respectively. Compared to the single-stage regression models, the validation errors of the two-stage regression models demonstrate no significant difference. For all fourteen VAVs in the validation set, the average RMSE and r2 of the predictions of the one-stage regression models are 15% more and 3% less than those of the two-stage regression models. These results show that the proposed models can indeed capture the underlying thermodynamics of the system.
3.6. Comparison with an ARMAX model According to Schunn and Wallach [49], one cannot apply absolute standards in assessing the quality of a particular goodness-offit value. Each goodness-of-fit should be compared to those obtained by previous models in the same domain. To this end, we shall compare the mpbARMAX1,2 models with existing models of the room temperature. Many linear parametric models such as ARX, ARMAX, BJ, OE, etc. have been successfully implemented to model the room temperature in HVAC units. In a recent work [24], an ARMAX model has been claimed to be very effective to predict the room temperature with relatively small errors in the summer and fall, the same periods when our data were collected. Following the modeling approach in [24], we test different combinations of na, nb, nc and nk, in which na, nb and nc range from 1 to 10 and nk ranges from 1 to 6. And we find that the ARMAX of order [na ¼ 1, nb ¼ 1, nc ¼ 1, nk ¼ 3], or ARMAX1113 for short, produces the best validation results on the SE building. This model is the same as one of the examples studied in [24]. The structure of the model ARMAX1113 is
Trm ðnÞ ¼
3.5. Model validation We first validate the pbARMAX1,2 models, which merely serves to assure that the appropriate models are chosen for the room temperature so that the mpbARMAX1,2 models are based on a solid foundation. We use the pbARMAX1,2 models to conduct four-week longterm predictions with the data from the training set. Take fourteen VAVs in the short wing on the third floor as examples. Tables 2and 3 in the previous Section 3.3 present the goodness-of-fits for the pbARMAX1,2 models. The average RMSE, MAE, r2, and MaxAE of the predictions are 0.12, 0.07, 0.97 and 1.26 for the pbARMAX1 model, while 0.10, 0.06, 0.98 and 1.24 for the pbARMAX2 model. Compared to the modeling errors presented in Section 3.3, the prediction demonstrates an increase in RMSE, and a slight decrease in r2. Nevertheless, the overall predictions of the pbARMAX1,2 models are accurate.
B1 ðqÞ B ðqÞ B ðqÞ _ nk Þ þ 2 Tdisch ðn nk Þ þ 3 Tclg ðn nk Þ mðn AðqÞ AðqÞ AðqÞ B4 ðqÞ B5 ðqÞ B ðqÞ þ Toa ðn nk Þ þ T ðn nk Þ þ 6 RHdisch AðqÞ AðqÞ htg AðqÞ B7 ðqÞ CðqÞ ðn nk Þ þ RHoa ðn nk Þ þ eðnÞ; ð15Þ AðqÞ AðqÞ
where Tclg and Thtg present the temperatures of cooling water, and heating water; RHdisch and RHoa denote the relative humidities of the discharge air and outside air; A(q), Bi(q) and C(q) are polynomials in q defined in Appendix. The variables in Eq. (15) that don’t appear in the mpbARMAX models are dropped in the comparison. Table 6 The four-week prediction errors generated by the mpbARMAX1 model for the example VAVs in the training set.
V351 V352 V353 V354
RMSE
MAE
r2
MaxAE
0.3541 0.1273 0.2591 0.4786
0.3204 0.0890 0.2330 0.4695
0.9752 0.9894 0.9888 0.9916
2.0186 1.2364 1.5983 1.3024
Room temperature (K)
S. Wu, J.Q. Sun / Building and Environment 56 (2012) 69e77
298
measurements training prediction
Training RMSE = 0.1098 r2 = 0.9303
297 296 Prediction RMSE = 0.1842 r2 = 0.9191
295 294 41
42
43
Time (days) Fig. 4. The short-term prediction of the two-stage regression mpbARMAX2 model over the validation set of VAVs demonstrates goodness-of-fits nearly as good as those over the training set.
Table 7 The four-week prediction errors generated by the mpbARMAX2 model for the example VAVs in the validation set.
V251 V252 V253 V254
RMSE
MAE
r2
MaxAE
0.3431 0.1462 0.2557 0.2839
0.3009 0.1035 0.1986 0.2755
0.9281 0.9800 0.9264 0.9930
1.4529 1.1708 2.1889 0.9272
The VAVs in the validation set are used for the comparison. The mpbARMAX1,2 models are formed as discussed previously. The ARMAX1113 model uses the historical data to train the model and obtains the optimal parameters. We apply the optimal training lengths for both approaches. We consider the prediction of the room temperature over a period of four weeks. Take V259 as an example. Fig. 5 compares the predictions from both models. The ARMAX1113 model performs well by providing a small RMSE as 0.0992, and a large r2 as 0.9712. The mpbARMAX1 model demonstrates an even better prediction with a 85.5% decrease in RMSE to 0.0143, and a 1.7% increase in r2 to 0.9878. Furthermore, the mpbARMAX1 model behaves better especially when the temperature is oscillating around the cooling setpoint. The observed MAE of the mpbARMAX1 model is 79.9% smaller than that of the ARMAX1113 model. For other thirteen VAVs in the validation set as well as the mpbARMAX2 model, we have observed the similar results. Table 8 compares the average goodness-of-fit of both the mpbARMAX1,2 model and the ARMAX1113 model for all VAVs in the
Fig. 5. The mpbARMAX1 model in the four-week long-term predictions over an example VAV outperforms the ARMAX1113 model in terms of all four validation metrics.
75
Table 8 The long-term prediction errors over four weeks by the mpbARMAX models (MR) and the ARMAX1113 model (AR). The mpbARMAX models demonstrate better performances consistently. MAE
r2
MaxAE
Prediction (weeks)
RMSE MR/AR
MR/AR
MR/AR
MR/AR
1 2 3 4
0.0420/0.0672 0.0541/0.1373 0.0591/0.1498 0.0908/0.1643 116%[/144%[
0.0382/0.0351 0.0384/0.0612 0.0431/0.0700 0.0806/0.0895 110%[/154%[
0.9995/0.9876 0.9980/0.9807 0.9974/0.9724 0.9959/0.9587 0.36%Y/3.92%Y
0.2414/0.6020 0.3335/1.1166 0.3984/1.1812 0.4853/1.7723 101%[/194%[
validation set. Compared to the ARMAX1113 model, the mpbARMAX1,2 models have a 19.5% smaller increase in RMSE over the four weeks, a 29.6% smaller increase in MAE, a 48.0% smaller increase in MaxAE, while a 90.9% less decrease in r2. This illustrates that the mpbARMAX1,2 models produce more accurate long-term predictions than the ARMAX1113 model. 4. Discussions and conclusions (1) Different training lengths in time for the first stage regression model lead to slightly different sets of coefficients a, b, c and d. The overall effect of such differences on the predictions of the first stage regression models is limited. Compared to the models identified with the optimal training lengths, the models identified with twenty days less or more than the optimal training length have a maximum 5% difference in goodness-of-fit of the prediction. (2) Our data were sampled at 5-min interval. All the results reported in this paper are based on the 5-min interval. Sampling time intervals, however, affect the outcome of the regression results. We have tested 1 min, 10 min, 15 min and 30 min in developing the mpbARMAX models. The shorter or longer sampling interval is achieved by linearly interpolating the original data or by deleting the extra data points. For sampling intervals longer than 5 min, the resulting mpbARMAX models have less than 5% change in the metrics for validation. The mpbARMAX models developed with the 5-min sampling interval on average produce the best predictions for all the VAVs. On the other hand, the mpbARMAX models developed with 1-min sampling interval produce the worse results compared to all other sampling intervals. (3) The proposed model is data-driven and can be updated when significant changes occur in HVAC control strategies, major equipment and operation profiles. One can even periodically update the model regardless, leading to an adaptive modeling strategy of HVAC systems. (4) The proposed model could also be used to predict other quantities of interest in HVAC systems such as the relative humidity. We limit the discussions to the temperature model in this paper only because the relative humidity measurements are unavailable. We are currently trying to manage a way to estimate the relative humidity based on the available measurements at the AHU level. We shall report this study in the future publications. To conclude, we have developed a new strategy to incorporate the architectural parameters into linear parametric regression models of the room temperature. The regression model is developed in two stages. In the first stage, linear parametric regression models of the room temperature based on the thermodynamic equations are developed with the data from VAVs involving different architectural parameters. The collection of the parameters
76
S. Wu, J.Q. Sun / Building and Environment 56 (2012) 69e77
of the regression models are then related to the architectural parameters with the help of another set of linear regression relationships. This is called the second stage regression. The model developed this way has been validated with the data in the validation set. We have found that the two-stage regression models are able to predict the room temperature in both short and long terms with a high accuracy. Furthermore, the proposed models demonstrate improved goodness-of-fit compared to a representative ARMAX model in the literature which has been proved to be very effective in the same application. Acknowledgment The authors would like to thank John Elliott of UC Merced for helping us to get the research started with access to the HVAC web control system. We deeply appreciate the generous technical guidance and input from John Elliott, Philip Haves and Michael Wetter of LBNL, and David Auslander of UC Berkeley. This paper was based on work of a CITRIS Seed Grant project. Appendix The parameters involved in the pbARMAX model are defined as follows.
! hwai Lwa k2wa 1þ ; P1 ¼ kwa þhwai Lwa kwa Lwa ðhwai þhwao ÞþL2wa hwai hwao P2 ¼
kwa hwao ; kwa ðhwai þ hwao Þ þ Lwa hwai hwao
hwdi Lwd P3 ¼ kwd þ hwdi Lwd
1þ
(16)
(17) !
k2wd
; kwd Lwd hwdi þ hwdo þ L2wd hwdi hwdo (18)
P4 ¼
P5 ¼
kwd hwdo ; kwd hwdi þ hwdo þ Lwd hwdi hwdo 1
ri VCv
(19)
;
(20)
_ p ; Prm ¼ P5 Swa hwai ðP1 1Þ þ hwdi zðP3 1Þ mC
(21)
_ p; Pdisch ¼ P5 mC
(22)
Poa ¼ P5 Swa hwai P2 þ hwdi zP4 ;
(23)
DE ¼ P5 f þ DEh;d :
(24)
The polynomials of the transfer functions in the ARMAX1113 model are given by
AðqÞ ¼ 1 þ a1 q1 þ a2 q2 þ . þ ana qna ; Bi ðqÞ ¼ 1 þ b1;i q1 þ b2;i q2 þ . þ bnb ;i qnb ;
(25) i ¼ 1; 2; .; 7 (26)
CðqÞ ¼ 1 þ c1 q1 þ c2 q2 þ . þ anc qnc ; where na, nb, and nc are the orders of the polynomials in q.
(27)
References [1] http://buildingsdatabook.eren.doe.gov/TableView.aspx?table¼1.1.1, U.S. residential and commercial buildings total primary energy consumption (quadrillion btu and percent of total). [2] http://buildingsdatabook.eren.doe.gov/TableView.aspx?table¼1.1.3, Buildings share of U.S. primary energy consumption (percent). [3] http://buildingsdatabook.eren.doe.gov/TableView.aspx?table¼1.1.9, Buildings share of U.S. electricity consumption (percent). [4] Schein J, Bushby ST, Castro NS, House JM. A rule-based fault detection method for air handling units. Energy Build 2006;38(12):1485e92. [5] Yoshida H, Kumar S. Rarx algorithm based model development and application to real time data for on-line fault detection in VAV AHU units. IBPSA Build Simulat 1999;99:161e8. [6] Fanger PO. Thermal comfort. Copenhagen: Danish Technical Press; 1970. [7] Fanger PO, Melikov AK, Hanzawa H, Ring J. Air turbulence and sensation of draught. Energy Build 1988;12(1):21e39. [8] American Society of Heating Refrigerating and Air Conditioning Engineers (ASHARE), Thermal environmental conditions for human occupancy (ASHARE, Standard 55-1992). [9] International Standards Organization (ISO), Moderate thermal environments: determination of the PMV and PPD indices and specification of the conditions for thermal comfort (ISO 7730:1994). [10] Ascione F, Bellia L, Capozzoli A, Minichiello F. Energy saving strategies in airconditioning for museums. Appl Thermal Eng 2009;29(4):676e86. [11] Balaras CA, Dascalaki E, Gaglia A. HVAC and indoor thermal conditions in hospital operating rooms. Energy Build 2007;39(4):454e70. [12] Yan D, Song F, Yang X, Jiang Y, Zhao B, Zhang X, et al. An integrated modeling tool for simultaneous analysis of thermal performance and indoor air quality in buildings. Build Environ 2008;43(3):287e93. [13] Wei S, Sun Y, Li M, Lin W, Zhao D, Shi Y, et al. Indoor thermal environment evaluations and parametric analyses in naturally ventilated buildings in dry season using a field survey and PMVe-PPDe model. Build Environ 2011;46(6): 1275e83. [14] Engdahl F, Johansson D. Optimal supply air temperature with respect to energy use in a variable air volume system. Energy Build 2004;36(3):205e18. [15] Yang IH, Kim KW. Prediction of the time of room air temperature descending for heating systems in buildings. Build Environ 2004;39(1):19e29. [16] Orosa JA. A new modelling methodology to control HVAC systems. Expert Syst Appl 2011;38(4):4505e13. [17] Tanimoto J, Hagishima A. State transition probability for the Markov model dealing with on/off cooling schedule in dwellings. Energy Build 2005;37(3): 181e7. [18] Pedrini A, Westphal FS, Lamberts R. A methodology for building energy modelling and calibration in warm climates. Build Environ 2002;37(8e9): 903e12. [19] Wetter M. Modelica-based modelling and simulation to support research and development in building energy and control systems. J Build Perform Simulation 2009;2(2):143e61. [20] Lowry G, Lee MW. Modelling the passive thermal response of a building using sparse BMS data. Appl Energy 2004;78(1):53e62. [21] Ríos-Moreno GJ, Trejo-Perea M, Castañeda-Miranda R, HernándezGuzmán VM, Herrera-Ruiz G. Modelling temperature in intelligent buildings by means of autoregressive models. Autom Construct 2007;16(5):713e22. [22] Peitsman HC, Bakker VE. Application of black-box models to HVAC systems for fault detection. ASHRAE Trans 1996;102(1):628e40. [23] Yiu JCM, Wang S. Multiple ARMAX modeling scheme for forecasting air conditioning system performance. Energy Convers Manag 2007;48(8): 2276e85. [24] Mustafaraj G, Chen J, Lowry G. Development of room temperature and relative humidity linear parametric models for an open office using BMS data. Energy Build 2010;42(3):348e56. [25] Déqué F, Ollivier F, Poblador A. Grey boxes used to represent buildings with a minimum number of geometric and thermal parameters. Energy Build 2000; 31(1):29e35. [26] Ghiaus C, Chicinas A, Inard C. Grey-box identification of air-handling unit elements. Control Eng Pract 2007;15:421e33. [27] Wen J, Smith TF. Development and validation of online models with parameter estimation for a building zone with VAV system. Energy Build 2007; 39(1):13e22. [28] Leephakpreeda T. Grey prediction on indoor comfort temperature for HVAC systems. Expert Syst Appl 2008;34(4):2284e9. [29] Tsilingiris PT. Wall heat loss from intermittently conditioned spaces e the dynamic influence of structural and operational parameters. Energy Build 2006;38(8):1022e31. [30] Devgan S, Jain AK, Bhattacharjee B. Predetermined overall thermal transfer value coefficients for composite, hot-dry and warm-humid climates. Energy Build 2010;42(10):1841e61. [31] Aste N, Angelotti A, Buzzetti M. The influence of the external walls thermal inertia on the energy performance of well insulated buildings. Energy Build 2009;41(11):1181e7. [32] Hwang RL, Shu SY. Building envelope regulations on thermal comfort in glass facade buildings and energy-saving potential for PMV-based comfort control. Build Environ 2011;46(4):824e34.
S. Wu, J.Q. Sun / Building and Environment 56 (2012) 69e77 [33] Sozer H. Improving energy efficiency through the design of the building envelope. Build Environ 2010;45(12):2581e93. [34] Korolija I, Marjanovic-Halburd L, Zhang Y, Hanby VI. Influence of building parameters and HVAC systems coupling on building energy performance. Energy Build 2011;43(6):1247e53. [35] Wang X, Eisenbrey J, Zeitz M, Sun JQ. Multi-stage regression analysis of acoustical properties of polyurethane foams. Sound Vib 2004;273(4e5): 1109e17. [36] Wang X, Sun JQ. On the fatigue analysis of non-Gaussian stress processes with asymmetric distribution. J Vib Acoust 2005;127(6):556e66. [37] Cunha JB, Couto C, Ruano AE. Real-time parameter estimation of dynamic temperature models for greenhouse environmental control. Control Eng Pract 1997;5(10):1473e81. [38] Jiménez MJ, Madsen H, Andersen KK. Identification of the main thermal characteristics of building components using MATLAB. Build Environ 2008; 43(2):170e80. [39] Bird RB, Stewart WE, Lightfoot EN. Transport phenomena. New York: Wiley; 2001. [40] Baierlein R. Thermal physics. Cambridge, England: Cambridge University Press; 2003. [41] Perrot P. A to Z of thermodynamics. Oxford, England: Oxford University Press; 1998.
77
[42] Baguley T. Beyond ANOVA: from repeated measures to multilevel models. In: Proceedings of PsyPAG mathematics statistics and computing workshop. Nottingham: Nottingham Trent University; 2008. [43] Ahn J. Beyond single equation regression analysis: path analysis and multistage regression analysis. Am J Pharm Educ 2002;66:37e42. [44] Kristjansson SD, Kircher JC, Webb AK. Multilevel models for repeated measures research designs in psychophysiology: an introduction to growth curve modeling. Psychophysiology 2007;44(5):728e36. [45] Bland JM, Altman DG. Statistics notes: correlations, regression, and repeated data. Br Med J 1994;308:896. [46] Mardia KV, Kent JT, Bibby JM. Multivariate analysis (probability and mathematical statistics). Waltham, Massachusetts: Academic Press; 1980. [47] Maxwell SE, Delaney HD. Designing experiments and analyzing data: a model comparison perspective. London: Routledge Academic; 2003. [48] Bishop CM. Neural networks for pattern recognition. Oxford: Clarendon Press; 1995. [49] Schunn CD, Wallach D. Evaluating goodness-of-fit in comparison of models to data. In: Tack W, editor. Psychologie der Kognition: Reden and vorträge anlässlich der emeritierung von Werner Tack. Saarbrueken, Germany: University of Saarland Press; 2005. p. 115e54. [50] Wu S, Sun JQ. Cross-level fault detection and diagnosis of building HVAC systems. Build Environ 2011;46(8):1558e66.