R. Kelley Pace

LREC Endowed Chair of Real Estate

E.J. Ourso College of Business Administration

Louisiana State University

Baton Rouge, LA 70803-6308

(225)-388-6238

FAX: (225)-334-1227

kelley@spatial-statistics.com

kelley@pace.am

Ronald Barry

Associate Professor of Statistics

Department of Mathematical Sciences

University of Alaska

Fairbanks, Alaska 99775-6660

(907)-474-7226

FAX: (907)-474-5394

FFRPB@aurora.alaska.edu

John M. Clapp

Professor of Finance

University of Connecticut

368 Fairfield Road, U-41RE

Storrs, CT 06269-2041

(860) 486-5057

FAX: (860) 486-0349

johnc@sbaserv.sba.uconn.edu

Mauricio Rodriguez

Assistant Professor of Finance

Texas Christian University

P.O. Box 32868

Fort Worth, TX 76129

(817) 921-7514

FAX: (817) 921-7227

mrodriguez@zeta.is.tcu.edu

This paper appeared in:

Abstract

Using 70,822 observations on housing prices during 1969-91 from Fairfax County Virginia, this manuscript demonstrates the substantial benefits obtained by modeling the spatial as well as the temporal dependence of the data. Specifically, the spatio-temporal autoregression with 12 variables reduced median absolute error by 37.35% relative to an indicator-based model with 26 variables. One-step ahead forecasts also document the improved performance of the proposed spatio-temporal model.

In addition, the manuscript illustrates techniques for rapidly computing the estimates and shows how to compute indices for any location.

1. Introduction

To say that the price of a parcel of real estate depends upon its location and recent market events seems all too obvious. However, the optimal way of incorporating spatial and temporal (spatio-temporal) dependencies into empirically feasible pricing models does not seem quite so obvious.

The use of indicator variables provides the easiest, although not necessarily optimal, way to model spatio-temporal dependencies. Due to the unequivocal ordering of time, one can easily determine the number of temporal indicator variables. For example, given data with a span of 20 years, one might have 20 annual indicator variables (given no intercept). Unfortunately, the lack of a unique, optimal arrangement (tessellation) of space impedes the formation of a simple set of spatial indicator variables. To capture the local structure of the market argues for employing relatively many indicator variables. To maintain parsimony argues for employing relatively few indicator variables. Modeling both space and time with interactions further raises the number of required indicators. For example, if 10 indicators for time and 20 indicators for space seems adequate for modeling each separately, the joint set would require 200 spatio-temporal indicator variables (given no intercept). Hence, it appears difficult to maintain parsimony and successfully capture local spatio-temporal effects.

To better capture the effect of both spatial and temporal information on real estate prices, while overcoming some of the problems associated with indicator variable models, we introduce a spatio-temporal model which uses information from nearby, recently sold properties in predicting the value of a given property. Instead of assuming that each region has its own effect modeled by a separate parameter as with the indicator variable based models, the spatio-temporal models assume that nearby properties ("comparables") have the same relation to the observations ("subjects") across the entire sample. These spatio-temporal models essentially generalize conventional hedonic regressions and the adjustment grid method (Pace and Gilley (forthcoming)). This approach leads naturally to a more parsimonious description of the data than the indicator variable approach.

This specific model falls under the very broad classification of spatio-temporal autoregressive (STAR) models. Applications to business and economic data are rare (c.f., Pfeifer and Bodily (1990)). However, these models have been applied to other areas such as imaging or hydrology. For example, Szummer and Picard (1996) use a STAR model to examine evolving images such as moving water or fire. Deutsch and Ramos (1986) have used them to predict river heights. See Cressie (1993) for more details on STAR models.

As a specific example of the two approaches, we estimated both an indicator variable model and a spatio-temporal model using data from Fairfax County, Virginia. The selected spatio-temporal model uses 12 variables while the indicator variable based model required 26 variables. Despite using fewer degrees-of-freedom, the spatio-temporal model showed strong improvements in goodness-of-fit relative to the indicator variable based model. For example, the median absolute errors went from 0.1478 under the indicator variable model to 0.0926 under the spatio-temporal model, a reduction of 37.35%.

In addition to strong statistical performance, spatio-temporal models provide other advantages. For example, since the estimated model provides a price surface that evolves over time, one can construct indices over time for any given location or location surfaces at any given point in time.

In what follows, section 2 introduces the design of the spatio-temporal model, section 3 applies the model to housing data from Fairfax County, Virginia, and section 4 concludes with the key results.

This section presents a simple spatio-temporal model. Section 2.1 develops the spatio-temporal model based upon a generalization of the traditional autoregressive error model, section 2.2 describes the spatio-temporal and temporal weight matrices, section 2.3 discusses estimation of the model, and section 2.4 delves into computational issues.

2.1 Spatio-Temporal Autoregressions

Assume the following autoregressive process,

where Y denotes the n by 1 vector of observations on the dependent variable, X denotes the n by k matrix of observations on the independent variables of interest, denotes a k by 1 vector of parameters, denotes an n by 1 vector of normal iid errors, and W denotes an n by n spatio-temporal weight matrix. The diagonal entries of W contain zeros to prevent each observation from predicting itself. Also, W contains only non-negative elements. In addition, we assume that each row in W sums to 1. Hence, W is a linear filter (Davidson and MacKinnon (1993, p. 691)).

For further structure, we assume the observations have been ordered according to time with the first row of W corresponding to the oldest observation. As Clapp, Dolde, and Tirtiroglu (1995) documented, a given transactions price impounds information (e.g., from changes in local taxes and public services) relevant to the pricing of neighboring properties. Thus, it seems eminently reasonable to assume that the sales price of a neighboring property will influence the subject property only if the neighboring sale is earlier in time. The time ordering of W (and its components) greatly simplifies the matrix multiplications, equation solutions, and determinant computations necessary to estimate the model.

Essentially, multiplying the random variables X, Y by filters these variables so the resulting transformed random variables , do not evince autocorrelation. In a temporal context, the operation reduces to taking the current values of Y and subtracting the previous values (or an average of the previous values) of Y scaled by a constant less than 1 (the autoregressive parameter). In a spatial context, the operation reduces to taking the values of Y at each location and subtracting the average of the surrounding values of Y scaled by a constant less than 1 (the spatial autoregressive parameter).

In spatio-temporal estimation, we wish to perform similar operations. Namely, we wish to take the current value of a variable at a location and subtract an average of past, surrounding values scaled by a constant less than 1. To implement this in a flexible way, we begin by partitioning W into a matrix S which specifies spatial relations among previous observations and T which specifies temporal relations among previous observations. Hence, the matrix T is a lag operator for regularly observed data. Also, S functions in space much like a lag operator does in time. The requirement of specifying relations only among previous observations and the ordering of the observations by time means W is a lower triangular matrix (barring ties in time).

We could implement the spatio-temporal filtering in a variety of different ways. For example, we could filter for space and time additively , where , represent autoregressive parameters related to space and time. We assume these parameters have absolute values of 1 or less. While conceptually simple, this tacitly assumes no interaction between the spatial and the temporal effects. Alternatively, we could filter for time first and then space or space first and then time . However, this presupposes we know which to filter first. On the other hand, we could generalize these by linearly combining both these forms. Equation represents a further generalization of the linear combination of the two forms.

The form in subsumes the linear combination of filtering for space first and then time or time first and then space. The matrices ST and TS will not generally have the same values. Hence, it seems reasonable to allow these variables their own autoregressive parameters, , . The flexibility of the filtering should allow to encompass the true model in many situations.

Substituting the form of the filter in into yields .

We could further generalize this and allow the parameters to vary for each form of lag. Hence, we write the generalized model as,

where Z denotes the n by p₁ matrix of observations on unlagged independent variables, X denotes the n by p₂ matrix of observations on the independent variables that will also have lags, denote five p₂ by 1 vectors of parameters, and denotes a p₁ by 1 vector of parameters on the unlagged variables. Essentially synthesizes the autoregressive distributed lag model of time series and the mixed regressive spatially autoregressive model of spatial econometrics (Ord (1975), Anselin (1988)).

Like these antecedent models, various restrictions yield different models of interest. For example, the restrictions yield the usual regression of Y on X while the restrictions , , , , and yield an autoregression in errors as previously described by . If one believes in the general to specific model strategy, one should estimate the general model and test these restrictions before adopting a specific alternative such as the autoregression in errors specification.

The above model explains the dependent variables as a function of the independent variables, the temporal effects which apply to all observations, and the effects which arise due to proximity in both time and space. Separating out the aggregate temporal and spatio-temporal effects allows for flexible specification of the global effects in time T versus the local effects in space and time in S. The spatial weight matrix, S, has a temporal element as well since it only has non-zero entries for previously sold observations. Incorporating the interaction terms ST and TS allows for the modeling of potentially compound spatio-temporal effects.

2.2. Specification of the Spatial and Temporal Weight Matrices

We can impose some additional structure on S and T. Housing data displays uneven density in both the temporal and spatial dimensions. A distance of one mile in the center of the city does not act the same as a distance of one mile in the suburbs. Similarly, transactions volume in housing markets can vary greatly over time. Real estate appraisers, who have great practical experience in predicting short-run housing prices, use a fixed number of neighbors in space ("comparables") and tend to pick a fixed number of neighbors in time as well. The use of ordinal distance and time tends to reduce the problems created by uneven data densities. For example, when the density of observations over time or space is very high, T or S averages over a short time interval or small spatial area. Similarly, when observations occur infrequently, T or S may average over a long time interval or large spatial area. Nearest neighbors have often been employed in density estimation as variable bandwidth smoothers (Silverman (1986)).

A heuristic approach can sometimes clarify the logistical complexities of spatio-temporal relations. Suppose there are six observations occurring at time periods 0,1,…,5 from a set of six locations a,b,…,f which lie on a line in that order. Further suppose =. The nearest neighbor to a is b, the nearest neighbors to c are a and b, and so forth. In this example, we assume one temporal lag and two spatial neighbors (). Suppose we have some temporally ordered observations on Y from different locations as presented below. Since we order over time, the observations in Y will not usually have a perfect spatial ordering. Because of the zeros in the first row, the observation at time 0 does not play a direct role but does help in forming the temporal and spatial lags. Since only previously transacted observations appear in S, T, these are lower triangular.

Obviously, STY and TSY do not have the same values. This means one will obtain different results by filtering for space first and time second or time first and space second. As there was a perfect linear relation between Y and location, SY for the last observations perfectly matches Y. It takes some observations before SY can perform adequately because at first only a few potential "comparables" exist. This leads to the decision of how many rows to discard.

Switching to implementation details, we used a geographic information system (GIS) to determine the decimal latitude and longitude of each house from its address, an operation known as geocoding. In the analysis, we used the unprojected decimal latitude and longitude coordinates to compute the Euclidean distance d_ij between every pair of observations j and i where (i>j). We subsequently sorted these distances and formed the set of individual neighbor matrices where represents the closest previously sold neighbor (shortest distance), represents the second previously sold neighbor (second shortest distance) and so on. The first rows of these matrices may have all zeros due to a lack of previously transacted neighbors. These very sparse matrices have a 1 in each row and contain 0s otherwise (apart from the initial rows). One can also restrict the past time interval used for neighbor selection. A priori, we decided that 15 neighbors (=15) should capture the vast majority of the spatial effects. Based upon both our priors and some very approximate preliminary fitting, we decided to restrict the construction of S to neighbors within five years in time. Hence, we do not use "comparables" which are more than five years old.

We computed the overall spatial matrix S via,

where weights the relative effect of the lth individual neighbor matrix. Hence, S depends upon the parameters and , as well as upon the underlying metric. Thus, imposes an autoregressive distributed lag structure on the spatial variables. By construction, each row in S sums to 1, S is lower triangular, and has zeros on the diagonal.

The use of the individual neighbor matrices, , greatly speeds up investigation of the sensitivity of the results to different forms of S. The individual neighbor matrices themselves require some expense in computation. However, reweighting these as in requires very little time.

Similarly, for each row of T, we give weight to the immediately prior observations. We enforce examination of the past by only looking at the lower triangle of T (recall the observations are sorted by time) which means observations can be nonzero only if i > j,

.

Note, the oldest observations correspond to the initial rows of T. Also, m_T can differ from m_S.

We form all the quantities (e.g., TY, SY, TX, SX, STX, TSX, STY, TSY) needed in computing the estimates for the entire sample. We subsequently drop some initial number of observations from these quantities for use in the actual estimation sample. We do this for two reasons. First, this allowed us to adjust m_T without confounding its effects by changing the size of the estimation sample (we examined the interval 0< m_T <1600). Second, without retaining some prior observations, the spatio-temporal estimator could perform poorly initially as it would have a very small selection of previously sold neighbors to use in the first predictions. For the data described later, we drop 1,600 of the initial observations which represents around three years of market data.

2.3 Estimation

For the model in , the profile log-likelihood function in the parameters becomes,

By the triangular nature of T and S and given these matrices have zeros on the diagonals, the matrices TS and ST are also triangular and have zeros on the diagonals. Hence, the determinant of (I-W) equals 1, and the log-determinant equals 0. Thus, maximizing the likelihood equates to minimizing the SSE via OLS. Lagged dependent variables can create bias problems for OLS in small samples. However, OLS consistently estimates the parameters in large samples (provided the errors are not autocorrelated).

2.4. Computational Considerations

Pace (1997) and Pace and Barry (1997a, b) have detailed some of the computational advantages of sparse data structures in spatial statistics. Sparse matrices have large proportions of zeros and special algorithms employ this fact to greatly accelerate computations and to save memory. The use of a fixed number of neighbors in time and space insures the sparsity of S and T. S will generally have a density of m_S/n, where m_S is the number of neighbors in space (e.g., closest 15 observations in space). Similarly, T will have average density m_T/n, where m_T is the number of neighbors in time. Hence, as n rises, S and T become progressively sparser. This greatly aids the computational feasibility of the spatio-temporal model.

The strictly lower triangular nature of T and S avoids the need to compute the determinant of a general matrix. This greatly reduces the complexity of calculating the spatio-temporal estimates.

As implemented, the spatio-temporal model possesses another major computational advantage over a straightforward execution of the stated model. For a large, active housing market T could have potentially thousands of non-zero entries in each row. This would make T rather cumbersome and potentially prohibitive to manipulate. In the spatio-temporal model, however, T does not appear by itself but only in combination with other variables (TY, TX, TSX, TSY). Thus, each column of these contains the running averages of the respective variable over time. Efficient linear filter routines exist for computing such quantities. For example, the Matlab function "filter" implements this operation.

The use of running sums and averages to handle large weighting matrices is analogous to the well-known Newey-West estimator that deals with heteroskedasticity and autocorrelation (Green (1993, p. 375-379)). The estimator uses a weighted average of products and cross-products from an OLS regression.

3. An Application to the Fairfax Housing Data

This section presents the application of the spatio-temporal modeling to housing data from Fairfax County, Virginia. The first part discusses the data, the second part details a traditional model of the data, the third part examines the general spatio-temporal model, the fourth part examines the performance of a more parsimonious model, the fifth part examines the predictive performance of the models, while the last part describes an index surface.

3.1. Fairfax Data

We began with 73,835 observations on housing transactions over 1966 through 1991 from Fairfax County Virginia. We initially discarded observations with (1) fewer than four rooms; (2) no baths; (3) more than 7 baths; (4) land area over 5 acres; (5) no land area; (6) more than 3 half-baths; (7) having a sales price of under $10,000 or over $1,000,000; and (8) sale dates outside of 1966-1992. This left 72,422 observations for the overall sample.

After computing Y, X, TY, TX, SY, STX, TSY, and so forth, we dropped the initial 1,600 observations from each of these quantities. This corresponds to around three years of data. We did this because of the lack of good neighbors for the initial observations. We used the remaining 70,822 observations as the actual estimation sample.

3.2 A Traditional Model

For comparison, we estimated a traditional OLS hedonic pricing model using indicator variables for the years 1969-1991. We used a double-log specification in Age, Land area Other Rooms, and Bathrooms, where Other Rooms equals the total number of rooms less Bathrooms. In all, this model employed 26 variables as shown in Table 1.

The model has an R² of 0.7927, all variables have the expected signs, and all seem quite significant. The time indicator variables show the expected pattern of rising values (except for the last two years). The log-likelihood was -305,552.13.

The results from the traditional model approximately match those in the literature. For example, in their study of the value of a view amenity, Rodriguez and Sirmans (1994) controlled for location by focusing their analysis to a small sub-area of Fairfax County. They report an adjusted R² of 0.729. In addition, Case, Pollakowski, and Wachter (1991) use the same data and one area in Fairfax County (Springfield) to investigate price indices. The largest reported R² on their data set of 14,617 observations equals 0.828, despite using 138 variables.

Even though the traditional model seems to perform similarly to other models in the literature, the residuals grossly violate independence as Figure 1 illustrates. The pairwise correlation in neighboring residuals varies between 0.44 and 0.25. This obvious dependence suggests the potential for using this information to improve estimation, a tact followed in the next section.

3.3. A General Spatio-temporal Model

We estimated the model in where Z contains two variables, an intercept and an index going from 1 to n and where X contains the variables log(Age), log(Land Area), log(Other Rooms), and log(Bathrooms) where Other Rooms equals Total Rooms less Bathrooms. The index should capture any deterministic trends relating to the sequence of the observations (the sequence is highly correlated with time). The general spatio-temporal model has 28 parameters (26 in the regression equation and the two non-linear parameters).

Table 2 shows the estimation results. The model has an R² of 0.8656, the variables log(Age), log(Land Area), log(Other Rooms), and log(Bathrooms) display the expected signs, and all of these as well as their spatially lagged versions are highly significant. Note, the extremely large t ratio of 174.52 for SY, the spatially lagged dependent variable. Based upon approximate fitting, we selected the values of 0.75 for and 650 for (corresponding roughly to two months in time).

The log-likelihood for the spatio-temporal model was -290,199.55, vastly larger than the one for the traditional model (a difference of 15,352.58).

Finally, the higher magnitude of coefficient estimate for the variable STY (-0.7902) relative to the variable TSY (-0.0324) suggests the need to filter first for time and subsequently for space and not vice versa for these data.

While one could examine the meaning and interpretation of the estimated coefficients in more depth, the number of such terms makes this somewhat difficult. Before proceeding further in this direction, it seems worthwhile to look at a simpler model to see what performance penalty exists from reducing the dimensionality of the model.

3.4. A Parsimonious Spatio-temporal Model

In searching for a parsimonious alternative model, the relatively low significance of the intercept and the index in the general model suggests that first differences over time might perform well. Table 3 presents the estimates for a model using first differences in the variables log(Age), log(Land Area), log(Other Rooms), and log(Bathrooms) along with their spatial lags coupled with the spatially lagged, temporally differenced dependent variable (I-T)Y. Based upon approximate fitting, we selected the values of 0.75 for and 650 for (corresponding roughly to two months in time).

The parsimonious model uses 12 variables (10 in the regression equation and two entering into the construction of S and T) as opposed to 28 for the general model. The resulting log-likelihood was -290,339.8. A statistically significant difference exists between the parsimonious and the general model based upon a likelihood ratio test. However, the difference seems small given the number of observations and the number of variables. Due to the relatively good performance of the parsimonious model, we adopted it for further examination.

3.5 Performance of the Models

Comparing the parsimonious spatio-temporal model to the traditional indicator-based model as in Figure 1, we see the spatio-temporal model residuals show low pairwise correlations among spatial neighbors. Essentially, it has used the spatial information to improve the estimates.

Table 4 shows the empirical residual (sample) statistics for the traditional model and for the parsimonious spatio-temporal model. The spatio-temporal model shows large improvements over the traditional approach. For example, median absolute errors went from 0.1478 to 0.0926, a 37.35% decrease.

We also computed the one-step ahead forecast errors (recursive residuals) for the parsimonious spatio-temporal model. As expected, the ex-sample spatio-temporal residuals are somewhat larger than the sample spatio-temporal residuals. However, the ex-sample residuals from the parsimonious model have a median absolute value of only 0.1014, an increase of 9.50% over the sample residuals but still 31.39% less than the traditional model sample residuals.

Figure 2 depicts these recursive residuals ordered over Time, Age, Land Area, Other Rooms, and Bathrooms. Under the null hypothesis of correct model specification, no discernible pattern should exist in any of these plots. In fact, the plots do not display any gross patterns. This adds to our confidence in the parsimonious spatio-temporal model.

3.6. An Index Surface

Due to the use of spatial and temporal lags, the estimated spatio-temporal models effectively create a spatial surface which evolves over time. Hence, at any point on the spatial surface we could separate out an index over time or for any given point in time we could separate out a spatial surface.

To examine the construction of a location specific index more closely, we write the equation for the estimated parsimonious spatio-temporal model at a given location.

At a given location, the spatial weight matrix averages the values at nearby locations for previously sold properties for any desired point in time, t. As time progresses, the entries in the rows of change as nearby properties ("comparables") sell or until these "comparables" become too old (over 60 months in the past). The matrices contains non-zero entries for properties recently sold prior to time t. Hence, provides the average price of houses immediately prior to time t. Similarly, provides the typical sold house’s characteristics immediately prior to time t. Finally, represents the archetypal house associated with the index (in this case , the average of all X).

This index has several interpretable components. The term represents the simple average of housing prices over the recent past. The term represents the value of the difference between the typical house’s characteristics and the recent average of house characteristics. The term focuses on value of the difference between the typical house’s characteristics and the recent average of house characteristics over the neighborhood. The term adjusts for the changes in price of recently sold, nearby properties.

Figure 3 displays potential time indices created by taking six arbitrary locations from the Fairfax area and shows how a particular property would fare over time. In addition, the overall market index appears as well (as given by the symbol ox). As expected, some locations showed relative increases, some showed stable paths, while others varied greatly over time.

Naturally, the precision of estimation for these paths depends upon the density of observations proximate in space and time. The estimation procedure optimized over the number of neighbors considered. Hence, smoother indices over space (more spatial neighbors) would reduce variance at the cost of location-specific bias. The use of too many neighbors could raise mean-squared error. As a stylized fact, better specification of the regression equation often reduces the spatial correlations among the errors. Hence, better specifications would probably require fewer neighbors and result in even more location specific indices.

4. Conclusion

Real estate prices evolve over both time and space. Traditional practice involves regression of sales price on property characteristics and indicator variables for time and location (e.g, zip code indicators). Unfortunately, this typically requires the inclusion of more variables in the models than desired on the basis of parsimony to yield residuals without visible temporal and spatial dependencies.

This paper follows a filtering or transformation approach to improve estimation. Specifically, one can filter the data for the temporal effects followed by filtering for the spatial effects, vice-versa, or combine the two approaches. Ideally, simple models applied to the filtered (transformed) data should not display the gross error dependencies found with the original data. The filtering approach leads to models with both temporal and spatial lags.

Applying the filtering approach to housing data from Fairfax Virginia resulted in large improvements in estimation. Specifically, the R² went from 0.7927 using the traditional model to 0.8656 using the spatio-temporal model, the log-likelihood rose by 15,352.58, and the median absolute errors fell by 37.35%. These gains persisted even when using one-step ahead forecast errors as opposed to sample errors. Other models in the literature, such as Model L in Case, Pollakowski, and Wachter (1991, p. 301) display poorer goodness-of-fit (R²=0.828) despite using more variables (k=138).

The results show that the log of sales price of any property is strongly influenced by the sales prices of previously sold, neighboring properties. While this is not surprising, the extent of the neighborhood dependence is very high in Fairfax County. In fact 80.82% of the temporal difference in log sales prices propagates from the nearest 15 neighbors (each weighted by a geometrically declining factor of 0.75) to the subject property. This suggest the short range spatial-temporal dependence is much greater than previously thought. Of course, the strong temporal dependence is entirely consistent with the large literature on ARIMA and related models.

The use of the spatial lags (based upon previous sales) means that the relative premium or discount at any given location relative to the overall market changes over time depending upon the sales of nearby properties. Hence, one can construct a temporal index for any given location. Also, one could construct a map of locational premia or discounts for any given time. In fact, the estimated model yields an evolving price surface.

Computationally, the spatio-temporal estimates use only the OLS estimator on suitably defined variables. Most of the effort comes from finding the spatio-temporal lags of the independent and dependent variables. However, the construction of these variables really just relies upon the criteria for comparable selection employed for many years by appraisers. The application of these techniques to this very large data set shows the feasibility such estimation — applications to smaller data sets should prove quite easy.

In conclusion, the large improvements in goodness-of-fit and the reduction in observed correlation among residuals of the spatio-temporal model relative to the traditional indicator-based model should earn it a place in the panoply of real estate statistical methods, especially since the computational difficulty of the spatio-temporal estimates does not appear large and the spatio-temporal model provides other benefits such as location specific indices.

Notes