The first and second moments of response variables for SDOF systems with hysteretic nonlinearity are obtained by a direct linearization procedure. This adaptation in the implementation of well-known statistical linearization methods, provides concise, model-independent linearization coefficients that are well-suited for numerical solution. The method may be applied to systems which incorporate any hysteresis model governed by a differential constitutive equation, and may be used for zero or non-zero mean random vibration. The implementation eliminates the effort of analytically deriving specific linearization coefficients for new hysteresis models. In doing so, the procedure of stochastic analysis is made independent from the task of physical modeling of hysteretic systems. In this study, systems with three different hysteresis models are analyzed under various zero and non-zero mean Gaussian White noise inputs. Results are shown to be in agreement with previous linearization studies and Monte Carlo Simulation.