Tutorial

Time Series Analysis & Forecasting in Python: Trend, Seasonality & Moving Averages

Complete guide to time series analysis and forecasting in Python. Master datetime indexing, resampling, moving average smoothing, seasonal decomposition, and the Augmented Dickey-Fuller (ADF) stationarity test.

Anuj SainiAug 24, 202612 min read

Every business decision involves anticipating what happens next: How much inventory should a warehouse stock next month? How many customer support agents should be scheduled for Friday? What will quarterly revenue look like?

Raw time series data is notoriously noisy. Spikes on weekends, holiday surges, and random daily fluctuations disguise long-term business trajectories.

In this guide, we'll build a complete Python time series analysis and forecasting framework from the ground up using Pandas, Seaborn, and Statsmodels.

For broader foundational Python data skills, check our Python Pandas Data Analysis Guide or follow the comprehensive Data Analyst Career Track.



Practice Time Series Aggregations

Master date parsing, interval windowing, and trailing averages on Topfolio Practice with instant browser-based execution.


1. The Anatomy of a Time Series

A time series is a sequence of observations recorded at successive, equally spaced points in time. Classically, it is composed of four components:

                           TIME SERIES DECOMPOSITION
               
     Observed Sales ($Y_t$)  =  Trend ($T_t$)  +  Seasonality ($S_t$)  +  Noise ($I_t$)
     
     1. TREND (T):       Long-term upward or downward trajectory over years.
     2. SEASONALITY (S): Fixed, repeating cyclical patterns (e.g., weekend spikes).
     3. NOISE (I):       Irregular, unpredictable stochastic fluctuations.
  • Additive Model: $Y_t = T_t + S_t + I_t$ (Used when seasonal variations are constant in size).
  • Multiplicative Model: $Y_t = T_t \times S_t \times I_t$ (Used when seasonal swings grow proportionally with sales volume).

2. Generating Realistic Time Series Data in Python

Let's simulate two full years of daily e-commerce sales with an upward linear trend, strong weekly seasonality, and Gaussian noise:

python
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
 
# Configure visual styling
sns.set_theme(style='whitegrid')
pd.set_option('display.float_format', lambda x: f'{x:.2f}')
 
# 1. Create a continuous daily timeline (2024 to 2025)
dates = pd.date_range(start='2024-01-01', end='2025-12-31', freq='D')
n = len(dates)
 
# 2. Additive Components
# Linear trend: Sales grow from $120/day to $380/day
trend = np.linspace(start=120, stop=380, num=n)
 
# Weekly seasonality: Weekend spikes (7-day periodic wave)
seasonality = 35 * np.sin(2 * np.pi * np.arange(n) / 7)
 
# Gaussian stochastic noise
np.random.seed(42)
noise = np.random.normal(loc=0, scale=18, size=n)
 
# Combine into observed revenue
sales = trend + seasonality + noise
 
# 3. Create DataFrame and enforce DatetimeIndex
df = pd.DataFrame({'Date': dates, 'Sales': sales})
df.set_index('Date', inplace=True)
 
print(f"Timeline: {df.index.min().strftime('%Y-%m-%d')} to {df.index.max().strftime('%Y-%m-%d')}")
print(f"Total Observations: {len(df):,} days")
df.head()
DateSales ($)
2024-01-01$128.94
2024-01-02$147.20
2024-01-03$153.11
2024-01-04$132.45
2024-01-05$112.80

2.1 Visualizing the Raw Time Series

python
plt.figure(figsize=(14, 5))
plt.plot(df.index, df['Sales'], color='#3498db', alpha=0.7, label='Observed Daily Sales')
plt.title("Daily E-Commerce Sales (Raw Noisy Series)", fontsize=14, fontweight='bold', pad=12)
plt.ylabel("Revenue ($)")
plt.xlabel("Date")
plt.legend(loc='upper left')
plt.tight_layout()
plt.show()

3. Resampling and Frequency Aggregations in Pandas

Daily data is often too volatile for executive reporting. Pandas provides the .resample() method to aggregate time series across different frequency cadences:

python
# Aggregate into Weekly Total and Monthly Mean
weekly_sales = df['Sales'].resample('W').sum()
monthly_avg_sales = df['Sales'].resample('MS').mean() # MS = Month Start
 
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(14, 4.5))
 
ax1.plot(weekly_sales.index, weekly_sales.values, color='purple', marker='o', markersize=3)
ax1.set_title("Weekly Total Sales Volume", fontweight='bold')
ax1.set_ylabel("Total Revenue ($)")
 
ax2.plot(monthly_avg_sales.index, monthly_avg_sales.values, color='darkgreen', marker='s', linewidth=2)
ax2.set_title("Monthly Average Daily Sales", fontweight='bold')
ax2.set_ylabel("Avg Daily Revenue ($)")
 
plt.tight_layout()
plt.show()
Frequency AliasMeaningTypical Usage
'D'Calendar DayDaily operational monitoring
'W' / 'W-MON'Weekly (Monday anchor)Sprint planning & weekly KPIs
'MS'Month StartFinancial accounting & cohort reporting
'QS'Quarter StartExecutive strategic forecasts

4. Signal from Noise: Moving Averages & Smoothing

A Moving Average (Rolling Mean) replaces each data point with the average of its surrounding window, dampening high-frequency noise to highlight underlying momentum.

4.1 Simple Moving Averages (7-Day and 30-Day)

python
# Calculate rolling 7-day and 30-day means
df['MA_7'] = df['Sales'].rolling(window=7, min_periods=1).mean()
df['MA_30'] = df['Sales'].rolling(window=30, min_periods=1).mean()
 
# Exponential Moving Average (14-day half-life span)
df['EMA_14'] = df['Sales'].ewm(span=14, adjust=False).mean()
 
plt.figure(figsize=(14, 6))
plt.plot(df.index, df['Sales'], label='Raw Daily Sales', color='gray', alpha=0.3)
plt.plot(df.index, df['MA_7'], label='7-Day Moving Avg (Weekly Trend)', color='#e67e22', linewidth=1.5)
plt.plot(df.index, df['MA_30'], label='30-Day Moving Avg (Macro Trend)', color='#27ae60', linewidth=2.5)
plt.title("Time Series Noise Filtering with Moving Averages", fontsize=14, fontweight='bold', pad=12)
plt.ylabel("Revenue ($)")
plt.legend(loc='upper left')
plt.tight_layout()
plt.show()

5. Classical Seasonal Decomposition with Statsmodels

To programmatically isolate the trend, seasonality, and residual noise, we use seasonal_decompose from statsmodels.tsa.seasonal:

python
from statsmodels.tsa.seasonal import seasonal_decompose
 
# Decompose with a 7-day seasonality period
decomposition = seasonal_decompose(df['Sales'], model='additive', period=7)
 
fig = decomposition.plot()
fig.set_size_inches(13, 8)
plt.suptitle("Additive Time Series Decomposition (Period = 7 Days)", fontsize=14, fontweight='bold', y=1.02)
plt.tight_layout()
plt.show()

5.1 Day-of-Week Seasonality Profiling

python
df['DayOfWeek'] = df.index.day_name()
day_order = ['Monday', 'Tuesday', 'Wednesday', 'Thursday', 'Friday', 'Saturday', 'Sunday']
 
avg_day_sales = df.groupby('DayOfWeek')['Sales'].mean().reindex(day_order)
 
plt.figure(figsize=(9, 4.5))
sns.barplot(x=avg_day_sales.index, y=avg_day_sales.values, palette='viridis')
plt.title("Average Sales by Day of Week (Weekly Seasonality)", fontsize=13, fontweight='bold', pad=10)
plt.ylabel("Avg Revenue ($)")
plt.xlabel("Day of Week")
plt.xticks(rotation=20)
plt.tight_layout()
plt.show()

6. Stationarity & The Augmented Dickey-Fuller (ADF) Test

A time series is stationary if its statistical properties (mean, variance, and autocorrelation) do not change over time. Most forecasting algorithms (ARIMA, SARIMA) strictly require stationary inputs.

6.1 Running the ADF Test

python
from statsmodels.tsa.stattools import adfuller
 
def test_stationarity(timeseries: pd.Series, title: str):
    print(f"\n--- Stationarity Test: {title} ---")
    result = adfuller(timeseries.dropna(), autolag='AIC')
    
    adf_stat = result[0]
    p_value = result[1]
    crit_values = result[4]
    
    print(f"ADF Statistic:       {adf_stat:.4f}")
    print(f"p-value:             {p_value:.6f}")
    print("Critical Values:")
    for key, val in crit_values.items():
        print(f"   {key}: {val:.3f}")
        
    if p_value < 0.05:
        print("✅ The series is stationary (Reject Null Hypothesis at 5% alpha).")
    else:
        print("❌ The series is NOT stationary (Fail to reject Null Hypothesis). Differencing required.")
 
# Test raw series
test_stationarity(df['Sales'], "Raw Sales Series")
--- Stationarity Test: Raw Sales Series ---
ADF Statistic:       -0.4215
p-value:             0.906478
Critical Values:
   1%: -3.439
   5%: -2.865
   10%: -2.569
❌ The series is NOT stationary (Fail to reject Null Hypothesis). Differencing required.

6.2 Achieving Stationarity via First-Order Differencing

Differencing subtracts the previous observation from the current observation: ΔY_t = Y_t - Y_(t-1).

python
df['Sales_Diff'] = df['Sales'].diff()
 
# Re-test differenced series
test_stationarity(df['Sales_Diff'], "First-Order Differenced Sales (1-Lag)")
--- Stationarity Test: First-Order Differenced Sales (1-Lag) ---
ADF Statistic:       -16.8452
p-value:             0.000000
✅ The series is stationary (Reject Null Hypothesis at 5% alpha).

7. Practical Baseline Forecasting Models

Before building complex deep learning networks (LSTM) or Prophet models, always establish a benchmark baseline.

7.1 Linear Trend-Slope Extrapolation (30-Day Forecast)

python
# Forecast 30 days beyond dataset end
last_date = df.index.max()
forecast_horizon = 30
 
# Estimate the recent 30-day slope of the 30-day Moving Average
recent_macro = df['MA_30'].dropna().tail(30)
slope = (recent_macro.iloc[-1] - recent_macro.iloc[0]) / 30
last_val = recent_macro.iloc[-1]
 
future_dates = pd.date_range(start=last_date + pd.Timedelta(days=1), periods=forecast_horizon, freq='D')
trend_forecast = [last_val + (slope * i) for i in range(1, forecast_horizon + 1)]
 
forecast_df = pd.DataFrame({'Date': future_dates, 'Forecast': trend_forecast}).set_index('Date')
forecast_df.head()

7.2 Holt-Winters Exponential Smoothing

For accurate forecasting with trend and seasonality, use Holt-Winters Exponential Smoothing from statsmodels:

python
from statsmodels.tsa.holtwinters import ExponentialSmoothing
 
# Train on historical data
train_model = ExponentialSmoothing(
    df['Sales'],
    trend='add',
    seasonal='add',
    seasonal_periods=7
).fit()
 
# Predict next 30 days
hw_forecast = train_model.forecast(steps=30)
forecast_df['Holt_Winters_Forecast'] = hw_forecast.values
 
# Plot the forecast
plt.figure(figsize=(14, 6))
plt.plot(df.index[-90:], df['Sales'][-90:], label='Actual Historical Sales (Last 90 Days)', color='#2c3e50', alpha=0.6)
plt.plot(df.index[-90:], df['MA_30'][-90:], label='30-Day Moving Average Trend', color='#27ae60', linestyle='--')
plt.plot(forecast_df.index, forecast_df['Holt_Winters_Forecast'], label='Holt-Winters 30-Day Forecast', color='#e74c3c', linewidth=2.5)
 
plt.title("Sales Forecast: 30-Day Future Projection", fontsize=14, fontweight='bold', pad=12)
plt.ylabel("Revenue ($)")
plt.legend(loc='upper left')
plt.tight_layout()
plt.show()

8. Forecast Accuracy Evaluation Metrics

To rigorously benchmark forecast models on historical test holds:

python
def calculate_forecast_metrics(actual: np.ndarray, predicted: np.ndarray) -> dict:
    mae = np.mean(np.abs(actual - predicted))
    rmse = np.sqrt(np.mean((actual - predicted) ** 2))
    mape = np.mean(np.abs((actual - predicted) / actual)) * 100
    
    return {
        'MAE ($)': mae,
        'RMSE ($)': rmse,
        'MAPE (%)': mape
    }
 
print("Forecast Evaluation Metrics Template Initialized.")

Complete End-to-End Time Series Script

You can run this complete diagnostic pipeline on any time-indexed CSV dataset:

python
import pandas as pd
from statsmodels.tsa.stattools import adfuller
from statsmodels.tsa.seasonal import seasonal_decompose
 
def analyze_time_series(df: pd.DataFrame, date_col: str, value_col: str, period: int = 7):
    """
    Automated time series inspection, smoothing, stationarity check, and decomposition.
    """
    df[date_col] = pd.to_datetime(df[date_col])
    ts = df.set_index(date_col)[value_col].asfreq('D').fillna(method='ffill')
    
    # 1. Moving Averages
    ma7 = ts.rolling(window=7).mean()
    ma30 = ts.rolling(window=30).mean()
    
    # 2. Stationarity
    adf_result = adfuller(ts.dropna())
    is_stationary = bool(adf_result[1] < 0.05)
    
    # 3. Decomposition
    decomp = seasonal_decompose(ts, model='additive', period=period)
    
    return {
        'series': ts,
        'ma7': ma7,
        'ma30': ma30,
        'adf_p_value': adf_result[1],
        'is_stationary': is_stationary,
        'decomposition': decomp
    }

Summary & What to Learn Next

Time series analysis allows data analysts to see past the noise of daily fluctuations and uncover the underlying drivers of business performance. By mastering datetime indices, resampling intervals, rolling window averages, seasonal decomposition, and the Augmented Dickey-Fuller stationarity test, you possess the essential toolkit to forecast demand and revenue.

Next Steps in Your Analytics Roadmap:

Frequently Asked Questions

What is the difference between Additive and Multiplicative time series decomposition?

In an Additive model (Y = Trend + Seasonality + Noise), seasonal fluctuations remain constant in magnitude regardless of the overall trend level. In a Multiplicative model (Y = Trend * Seasonality * Noise), the amplitude of seasonal swings expands proportionally as the trend grows (common in retail sales and revenue growth).

Why must a time series be stationary before fitting statistical models like ARIMA?

Statistical forecasting models require stationarity—constant mean, constant variance, and autocovariance independent of time—so that historical statistical properties remain valid for future prediction. Non-stationary series with trends or variable variance produce spurious regressions unless transformed via differencing or log scaling.

How does .resample() differ from .rolling() in Pandas?

.resample() changes the temporal granularity of the dataset by grouping and aggregating timestamps into fixed buckets (e.g., converting 365 daily rows into 12 monthly totals). .rolling() maintains the original row count and computes sliding window summary statistics (e.g., 7-day trailing average) across a moving window.

What does a p-value < 0.05 in the Augmented Dickey-Fuller (ADF) test indicate?

The null hypothesis (H0) of the ADF test states that a unit root is present (the time series is non-stationary). A p-value < 0.05 rejects H0 at the 5% significance level, providing strong statistical evidence that the series is stationary and ready for autoregressive modeling.

What evaluation metrics should I use to compare time series forecast models?

Use Mean Absolute Error (MAE) for intuitive average error in original units, Root Mean Squared Error (RMSE) when large prediction errors carry severe penalties, and Mean Absolute Percentage Error (MAPE) to express error as a relative percentage of actual values across different scales.

Anuj Saini

Written by

Anuj SainiFounder & Lead Instructor

Founder at Topfolio with 6+ years in data & analytics across JPMC, Ultrahuman, and high-growth startups. Sat on hiring panels, reviewed 500+ resumes, and writes practical SQL & data guides.