A One-Dimensional (1D) Reduced-Order Model (ROM) has been developed for a 3D Rayleigh-B\'enard convection system in the turbulent regime with Rayleigh number $\mathrm{Ra}=10^6$. The state vector of the 1D ROM is horizontally averaged temperature. Using the Green's Function (GRF) method, which involves applying many localized, weak forcings to the system one at a time and calculating the responses using long-time averaged Direct Numerical Simulations (DNS), the system's Linear Response Function (LRF) has been computed. Another matrix, called the Eddy Flux Matrix (EFM), that relates changes in the divergence of vertical eddy heat fluxes to changes in the state vector, has also been calculated. Using various tests, it is shown that the LRF and EFM can accurately predict the time-mean responses of temperature and eddy heat flux to external forcings, and that the LRF can well predict the forcing needed to change the mean flow in a specified way (inverse problem). The non-normality of the LRF is discussed and its eigen/singular vectors are compared with the leading Proper Orthogonal Decomposition (POD) modes of the DNS data. Furthermore, it is shown that if the LRF and EFM are simply scaled by the square-root of Rayleigh number, they perform equally well for flows at other $\mathrm{Ra}$, at least in the investigated range of $ 5 \times 10^5 \le \mathrm{Ra} \le 1.25 \times 10^6$. The GRF method can be applied to develop 1D or 3D ROMs for any turbulent flow, and the calculated LRF and EFM can help with better analyzing and controlling the nonlinear system.