This paper considers panel data models with cross-sectional dependence arising from both spatial autocorrelation and unobserved common factors. It derives conditions for model identification and proposes estimation methods that employ cross-sectional averages as factor proxies, including the 2SLS, Best 2SLS, and GMM estimations. The proposed estimators are robust to unknown heteroskedasticity and serial correlation in the disturbances, unrequired to estimate the number of unknown factors, and computationally tractable. The paper establishes the asymptotic distributions of these estimators and compares their consistency and efficiency properties. Extensive Monte Carlo experiments lend support to the theoretical findings and demonstrate the satisfactory finite sample performance of the proposed estimators. The empirical section of the paper finds strong evidence of spatial dependence of real house price changes across 377 Metropolitan Statistical Areas in the US from 1975Q1 to 2014Q4. The results also reveal that population and income growth have significantly positive direct and spillover effects on house price changes. These findings are robust to different specifications of the spatial weights matrix constructed based on distance, migration flows, and pairwise correlations. . This work was carried out during my doctoral study at the University of Southern California.(ii) The error term e it has absolutely summable cumulants up to the fourth order.Assumption 3. The factor loadings, γ i and A i , are independently and identically distributed across i, and independent of e jt , v jt , and f t , for all i, j, and t. Both γ i and A i have fixed means, which are given by γ and A, respectively, and finite variances. In particular, for all i,where Ω η is a symmetric non-negative definite matrix, γ < K, A < K, andAssumption 4. The true parameter vector, δ 0 = ρ 0 , β ′ 0 ′ , is in the interior of the parameter space, denoted by ∆ sp , which is a compact subset of the (k + 1)-dimensional Euclidean space, R k+1 .Assumption 5. The matrixC, given by (7), has full row rank for all N , including N → ∞. Assumption 6. The N × N nonstochastic spatial weights matrix, W = (w ij ), has bounded row and column sum norms, namely, ||W|| ∞ < K and ||W|| 1 < K, respectively, andfor all values of ρ. In addition, the diagonal entries of W are zero, that is, w ii = 0, for all i = 1, 2, . . . , N . Assumption 7. The N × q matrix of instrumental variables, Q .t , for t = 1, 2, . . . , T , is composed of a subset of the columns of X .t , WX .t , W 2 X .t , . . . , and its column dimension q is fixed for all N and t. The matrix Q = (Q ′ .1 , Q ′ .2 , . . . , Q ′ .T ) ′ represents the IV matrix of dimension N T × q. Assumption 8. (i) There exists N 0 and T 0 , such that for all N > N 0 and T > T 0 , the matrices (N T ) −1 Q ′ M b Q and (N T ) −1 Q ′ M b f Q exist and are nonsingular.