Summary This paper presents two methods for generation artificial acceleration time series focused on the simulation of non‐Gaussian spiky fluctuation features. The random acceleration time series were recorded in railway boxcar; the simulations were performed under the PSD and higher statistical moments (skewness and kurtosis) restrictions. The first method is based on discrete Fourier transform method with EARPG(1) model, which, in nature, is a Fourier representation of time series with special selection of random Fourier phases. The synthetic time series generation are presented for a given target PSD property, skewness, and kurtosis. In the second method, we characterized the stationary and non‐Gaussian process using polynomial chaos expansion and Karhunen‐Loeve expansion approach based on the target PSD, skewness, and kurtosis. The accuracy of the methods are verified by moments and PSD comparison. The package is simplified to be a SDOF system; the acceleration response characteristics of the system are analyzed; the results can be applied in packaging system behaviors evaluation and packaging design optimization.