在R语言中实现时间序列平稳化并构建ARMA模型(新手求助)
Hey there! Since you're new to R and time series analysis, let's walk through this process step by step—from loading your raw data all the way to building and validating an ARMA model. I'll keep things simple with actionable code snippets you can copy and tweak for your data.
First, we'll turn your raw numerical data into a proper time series object in R. I'll assume your data is monthly (adjust the frequency argument if it's quarterly, yearly, etc.):
# Paste your full data into this vector (I included the snippet you provided) data_vector <- c(2953, 2635, 2404, 2413, 2136, 1565, 1451, 2037, 2477, 2785, 2994, 2681, 3098, 2708, 2517, 2445, 2087, 1801, 1216, 2173, 2286, 3121, 3458, 3511, 3524, 2767, 2744, 2603, 2527, 1846, 1066, 2327, 3066, 3048, 3806, 4042, 3583, 3438, 2957, 2885, 2744, 1837, 1447, 2504, 3248, 3098, 4318, 3561, 3316, 3379, 2717, 2354, 2445, 1542, 1606, 2590, 3588, 3202, 4704, 4005, 3810, 3488, 2781, 2944, 2817, 1960, 1937, 2903, 3357, 3552, 4581, 3905, 4581, 4037, 3345, 3175, 2808, 2050, 1719, 3143, 3756, 4776, 4540, 4309, 4563) # Convert to a time series object (frequency = 12 for monthly data) ts_data <- ts(data_vector, frequency = 12)
ARMA models require a stationary series (constant mean, variance, and autocorrelation over time). We'll use two tools to check this:
- A visual plot of the series
- The Augmented Dickey-Fuller (ADF) statistical test
First, install and load the necessary packages if you haven't already:
install.packages(c("tseries", "forecast")) library(tseries) library(forecast)
Now run the checks:
# Plot the raw series to spot trends/seasonality plot(ts_data, main = "Raw Time Series", ylab = "Value") # ADF test: Null hypothesis = series is non-stationary adf_result <- adf.test(ts_data) print(adf_result) # Plot ACF/PACF to see autocorrelation patterns par(mfrow = c(1, 2)) acf(ts_data, main = "ACF of Raw Series") pacf(ts_data, main = "PACF of Raw Series") par(mfrow = c(1, 1))
What to look for:
- If the ADF test p-value is less than 0.05, we reject the null hypothesis—your series is stationary.
- If the ACF plot shows slow decay (instead of dropping off quickly), that's a sign of non-stationarity.
If your series isn't stationary, we'll fix it with common techniques:
Option 1: First/Second Differencing
Differencing removes trends by calculating the difference between consecutive values:
# First difference (subtract each value from the next) diff1 <- diff(ts_data, differences = 1) plot(diff1, main = "First-Differenced Series") # Re-test stationarity adf_diff1 <- adf.test(diff1) print(adf_diff1) # Check ACF/PACF of differenced data par(mfrow = c(1, 2)) acf(diff1, main = "ACF of First-Differenced Series") pacf(diff1, main = "PACF of First-Differenced Series") par(mfrow = c(1, 1))
If first differencing isn't enough, try second differencing with differences = 2.
Option 2: Log Transformation + Differencing
If your series has increasing variance, take the log first before differencing (only works if all values are positive, which yours are):
log_ts <- log(ts_data) diff_log <- diff(log_ts, differences = 1) plot(diff_log, main = "Log + First-Differenced Series") adf_diff_log <- adf.test(diff_log) print(adf_diff_log)
Repeat until your ADF p-value is < 0.05.
ARMA models have two parameters:
p: Number of autoregressive termsq: Number of moving average terms
Manual Selection (using ACF/PACF):
- If the ACF cuts off sharply at lag q and PACF decays slowly → use MA(q)
- If the PACF cuts off sharply at lag p and ACF decays slowly → use AR(p)
- If both decay slowly → use ARMA(p,q)
Automatic Selection (easier for beginners):
Use auto.arima() to let R pick the best model for you. It will even handle differencing automatically (so you can pass the original series):
# auto.arima will test different p/d/q combinations (d = differencing order) best_model <- auto.arima(ts_data, seasonal = FALSE, trace = TRUE) print(best_model)
seasonal = FALSEassumes no seasonal pattern (set toTRUEif you see seasonality in your plot)trace = TRUEshows the model search process so you can follow along
If you chose orders manually, use the arima() function (set d=0 since we already stationarized the series):
# Example: Fit ARMA(1,1) (replace p=1 and q=1 with your chosen values) arma_model <- arima(ts_data, order = c(1, 0, 1)) print(arma_model)
If you used auto.arima(), the output already includes your fitted model.
A good ARMA model should have white noise residuals (no remaining autocorrelation). Let's check:
# Plot residual diagnostics (look for random scatter in residuals, no patterns) checkresiduals(best_model) # Ljung-Box test: Null hypothesis = residuals are white noise ljung_box <- Box.test(best_model$residuals, type = "Ljung-Box") print(ljung_box)
If the Ljung-Box p-value is greater than 0.05, we can't reject the null hypothesis—your residuals are white noise, meaning the model captures all predictable patterns in the data.
Once you have a valid model, you can forecast future values:
# Forecast next 12 periods (adjust h to your desired number of forecasts) forecasts <- forecast(best_model, h = 12) plot(forecasts)
内容的提问来源于stack exchange,提问作者makome

