Time Series Analysis

Time series analysis involves studying data collected over time to identify patterns, trends, seasonality, and anomalies. In Earth observation, this approach is powerful for monitoring environmental changes, agricultural cycles, and long-term trends.

The first step in time series analysis is collecting consistent data over time:

Step 1: Search for Available Imagery

import requests
import json
import pandas as pd
import matplotlib.pyplot as plt
from datetime import datetime

api_key = "your_api_key_here"
url = "https://observearth.com/api/s2/search/"

payload = {
    "geometry_id": "123e4567-e89b-12d3-a456-426614174000",
    "start_date": "2020-01-01",
    "end_date": "2022-12-31",
    "cloud_cover": 20
}

headers = {
    "X-API-Key": api_key,
    "Content-Type": "application/json"
}

response = requests.post(url, headers=headers, data=json.dumps(payload))
data = response.json()

print(f"Found {data['count']} images")

Step 2: Calculate Statistics for Each Image

url = "https://observearth.com/api/s2/stats/"

payload = {
    "geometry_id": "123e4567-e89b-12d3-a456-426614174000",
    "start_date": "2020-01-01",
    "end_date": "2022-12-31",
    "cloud_cover": 20,
    "index": "ndvi"
}

response = requests.post(url, headers=headers, data=json.dumps(payload))
stats_data = response.json()

# Create a time series dataframe
dates = []
ndvi_means = []

for result in stats_data['results']:
    dates.append(datetime.strptime(result['date'], '%Y-%m-%d'))
    ndvi_means.append(result['mean_value'])

time_series_df = pd.DataFrame({
    'date': dates,
    'ndvi_mean': ndvi_means
})

time_series_df = time_series_df.sort_values('date')
time_series_df.set_index('date', inplace=True)

Visualization is key to understanding patterns in time series data:

Basic Time Series Plot

# Plot the time series
plt.figure(figsize=(12, 6))
plt.plot(time_series_df.index, time_series_df['ndvi_mean'], marker='o', linestyle='-', markersize=4)
plt.title('NDVI Time Series (2020-2022)')
plt.xlabel('Date')
plt.ylabel('Mean NDVI')
plt.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()

Advanced Visualizations

# Resample to monthly data and plot with seasonal patterns
monthly_data = time_series_df.resample('M').mean()

# Create a figure with multiple years overlaid
years = monthly_data.index.year.unique()
plt.figure(figsize=(12, 6))

for year in years:
    year_data = monthly_data[monthly_data.index.year == year]
    plt.plot(year_data.index.month, year_data['ndvi_mean'], marker='o', linestyle='-', label=str(year))

plt.title('Monthly NDVI Patterns by Year')
plt.xlabel('Month')
plt.ylabel('Mean NDVI')
plt.xticks(range(1, 13), ['Jan', 'Feb', 'Mar', 'Apr', 'May', 'Jun', 'Jul', 'Aug', 'Sep', 'Oct', 'Nov', 'Dec'])
plt.legend(title='Year')
plt.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()

One of the most valuable aspects of time series analysis is the ability to detect long-term trends and identify anomalies:

Trend Analysis

from scipy import stats

# Calculate trend (linear regression)
x = np.arange(len(monthly_data))
y = monthly_data['ndvi_mean'].values
slope, intercept, r_value, p_value, std_err = stats.linregress(x, y)

# Create trend line
trend_line = intercept + slope * x

# Plot data with trend line
plt.figure(figsize=(12, 6))
plt.plot(monthly_data.index, monthly_data['ndvi_mean'], marker='o', linestyle='-', alpha=0.7, label='NDVI')
plt.plot(monthly_data.index, trend_line, 'r--', linewidth=2, label=f'Trend (slope={slope:.6f})')
plt.title('NDVI Trend Analysis')
plt.xlabel('Date')
plt.ylabel('Mean NDVI')
plt.legend()
plt.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()

print(f"Trend slope: {slope:.6f} per month")
print(f"P-value: {p_value:.6f}")
print(f"R-squared: {r_value**2:.6f}")

Anomaly Detection

# Calculate rolling mean and standard deviation
window = 12  # 12-month window
rolling_mean = monthly_data['ndvi_mean'].rolling(window=window, center=True).mean()
rolling_std = monthly_data['ndvi_mean'].rolling(window=window, center=True).std()

# Define anomalies (values outside 2 standard deviations)
threshold = 2
upper_bound = rolling_mean + (threshold * rolling_std)
lower_bound = rolling_mean - (threshold * rolling_std)

# Identify anomalies
anomalies = monthly_data[(monthly_data['ndvi_mean'] > upper_bound) | 
                         (monthly_data['ndvi_mean'] < lower_bound)].copy()

# Plot with anomalies highlighted
plt.figure(figsize=(12, 6))
plt.plot(monthly_data.index, monthly_data['ndvi_mean'], label='NDVI')
plt.plot(monthly_data.index, rolling_mean, 'r--', label='Rolling Mean')
plt.plot(monthly_data.index, upper_bound, 'g-.', alpha=0.5, label='Threshold (±2σ)')
plt.plot(monthly_data.index, lower_bound, 'g-.', alpha=0.5)
plt.scatter(anomalies.index, anomalies['ndvi_mean'], color='red', s=50, label='Anomalies')
plt.title('NDVI Anomaly Detection')
plt.xlabel('Date')
plt.ylabel('Mean NDVI')
plt.legend()
plt.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()

print(f"Detected {len(anomalies)} anomalies")

Change detection identifies significant changes between different time periods:

Change Detection Methods:

• Image differencing: Subtract one time period from another
• Ratio analysis: Calculate the ratio between two time periods
• Regression analysis: Model the relationship between time periods
• Classification comparison: Compare land cover classifications
• Change vector analysis: Analyze changes in multiple dimensions

Now that you understand time series analysis, you can:

• Learn about batch processing for large-scale analysis
• Explore integration with GIS and data science tools
• Work with Landsat's long-term archive for historical analysis