This paper presents a mathematically rigorous, subspace projection-based reduced-order modeling methodology and an integrated framework to automatically generate reduced-order models for spacecraft thermal analysis. Two key steps in the reduced-order modeling procedure are described: first, the acquisition of a full-scale spacecraft model in the ordinary differential equation, and differential algebraic equation, forms to resolve its dynamic thermal behavior; and second, reduced-order modeling to markedly reduce the dimension of the full-scale model. Specifically, proper orthogonal decomposition in conjunction with a discrete empirical interpolation method and trajectory piecewise-linear methods are developed to address the strong nonlinear thermal effects due to coupled conductive and radiative heat transfer in the spacecraft environment. Case studies using NASA-relevant satellite models are undertaken to verify the capability and to assess the computational performance of the reduced-order modeling technique in terms of speedup and error relative to the full-scale model. Reduced-order modeling exhibits excellent agreement in spatiotemporal thermal profiles (less than 0.5% relative error in pertinent timescales) along with salient computational acceleration (up to two orders of magnitude speedup) over the full-scale analysis. These findings establish the feasibility of reduced-order modeling to perform rational and computationally affordable thermal analysis, develop reliable thermal control strategies for spacecraft, and greatly reduce the development cycle times and costs. Nomenclature A = conduction exchange matrices C = thermal capacitance matrix D = scatter matrix of thermal loads Q = heat flux R = radiation exchange matrices T = temperature vector, K t = time, s U = projection space u = independent thermal inputs Subscripts r = quantity in reduced dimension space s = thermal links between the nodes and the isothermal objects t = time-varying quantity