Sample Data

def simple_forecast() -> xr.DataArray:
    """
    Generate simple series of prediction values for tutorials.

    Returns:
        DataArray with simple forecast values

    Characteristics:
        - Single dimension array
        - Realistic forecast values
        - No missing data
        - Suitable for basic scoring function demos
    """

def simple_observations() -> xr.DataArray:
    """
    Generate simple series of observation values for tutorials.

    Returns:
        DataArray with simple observation values

    Characteristics:
        - Matches simple_forecast() structure
        - Corresponding observation values
        - Suitable for basic verification examples
    """

def simple_forecast_pandas() -> pd.Series:
    """
    Generate simple pandas series of prediction values.

    Returns:
        Pandas Series with forecast values

    Notes:
        - Pandas Series format instead of xarray
        - Compatible with pandas-specific scoring functions
        - Useful for traditional pandas workflows
    """

def simple_observations_pandas() -> pd.Series:
    """
    Generate simple pandas series of observation values.

    Returns:
        Pandas Series with observation values

    Notes:
        - Matches simple_forecast_pandas() structure
        - Corresponding observation data
        - Suitable for pandas-based verification
    """

def continuous_observations(*, large_size: bool = False) -> xr.DataArray:
    """
    Create continuous observation array with synthetic data.

    Args:
        large_size: Generate larger dataset for performance testing

    Returns:
        Multi-dimensional DataArray with synthetic observation data

    Dimensions:
        - time: Temporal dimension with regular intervals
        - station: Spatial stations (if multi-dimensional)
        - Additional dimensions based on configuration

    Characteristics:
        - Realistic temporal and spatial patterns
        - Seasonal cycles and trends
        - Missing data patterns
        - Labeled coordinates with metadata
    """

def continuous_forecast(*, large_size: bool = False, lead_days: bool = False) -> xr.DataArray:
    """
    Create continuous forecast array with synthetic data.

    Args:
        large_size: Generate larger dataset for performance testing
        lead_days: Include lead time dimension for forecast horizons

    Returns:
        Multi-dimensional DataArray with synthetic forecast data

    Dimensions:
        - time: Valid time dimension
        - station: Spatial stations (if multi-dimensional)
        - lead_time: Forecast lead times (if lead_days=True)

    Characteristics:
        - Corresponds to continuous_observations() structure
        - Realistic forecast errors and biases
        - Lead time dependencies (if enabled)
        - Ensemble-like variations
    """

def cdf_forecast(*, lead_days: bool = False) -> xr.DataArray:
    """
    Create forecast array with CDF at each point.

    Args:
        lead_days: Include lead time dimension

    Returns:
        DataArray with CDF forecast values [0, 1]

    Dimensions:
        - time: Valid time dimension
        - threshold: CDF threshold values
        - lead_time: Forecast lead times (if lead_days=True)

    Characteristics:
        - Monotonically increasing CDFs
        - Realistic probability distributions
        - Multiple threshold levels
        - Suitable for CRPS calculations
    """

def cdf_observations() -> xr.DataArray:
    """
    Create observation array compatible with CDF forecasts.

    Returns:
        DataArray with observation values

    Characteristics:
        - Compatible with cdf_forecast() output
        - Continuous values for CDF evaluation
        - Matching temporal structure
        - Realistic value ranges
    """

from scores.sample_data import simple_forecast, simple_observations
from scores.continuous import mse, rmse, mae

# Generate simple tutorial data
forecast = simple_forecast()
observations = simple_observations()

print(f"Forecast shape: {forecast.shape}")
print(f"Forecast range: [{forecast.min().values:.2f}, {forecast.max().values:.2f}]")

# Calculate basic scores
mse_score = mse(forecast, observations)
rmse_score = rmse(forecast, observations)
mae_score = mae(forecast, observations)

print(f"\nBasic Scores:")
print(f"MSE:  {mse_score.values:.3f}")
print(f"RMSE: {rmse_score.values:.3f}")
print(f"MAE:  {mae_score.values:.3f}")

from scores.sample_data import simple_forecast_pandas, simple_observations_pandas
from scores.pandas import mse as pandas_mse

# Generate pandas data
forecast_pd = simple_forecast_pandas()
observations_pd = simple_observations_pandas()

print(f"Pandas forecast type: {type(forecast_pd)}")
print(f"Pandas data length: {len(forecast_pd)}")

# Use pandas-specific scoring functions
mse_pd = pandas_mse(forecast_pd, observations_pd)
print(f"Pandas MSE: {mse_pd:.3f}")

from scores.sample_data import continuous_forecast, continuous_observations
from scores.continuous import mse, kge, pearsonr

# Generate multi-dimensional data
forecast = continuous_forecast()
observations = continuous_observations()

print(f"Multi-dimensional data:")
print(f"Forecast dimensions: {forecast.dims}")
print(f"Forecast shape: {forecast.shape}")
print(f"Coordinates: {list(forecast.coords.keys())}")

# Analyze different aspects
temporal_mse = mse(forecast, observations, reduce_dims="time")
spatial_mse = mse(forecast, observations, preserve_dims="time")
overall_mse = mse(forecast, observations)

print(f"\nMulti-dimensional Analysis:")
print(f"Temporal MSE shape: {temporal_mse.shape}")
print(f"Spatial MSE shape: {spatial_mse.shape}")
print(f"Overall MSE: {overall_mse.values:.3f}")

# Advanced metrics
kge_score = kge(forecast, observations)
correlation = pearsonr(forecast, observations)

print(f"KGE: {kge_score.values:.3f}")
print(f"Correlation: {correlation.values:.3f}")

# Generate large datasets for performance testing
large_forecast = continuous_forecast(large_size=True)
large_observations = continuous_observations(large_size=True)

print(f"Large dataset dimensions: {large_forecast.shape}")
print(f"Memory usage estimate: ~{large_forecast.nbytes / 1e6:.1f} MB")

# Performance timing example
import time

start_time = time.time()
large_mse = mse(large_forecast, large_observations)
end_time = time.time()

print(f"Large dataset MSE: {large_mse.values:.3f}")
print(f"Computation time: {end_time - start_time:.3f} seconds")

from scores.sample_data import continuous_forecast

# Generate forecast with lead times
lead_forecast = continuous_forecast(lead_days=True)

print(f"Lead time forecast dimensions: {lead_forecast.dims}")
print(f"Lead time coordinate: {lead_forecast.lead_time.values}")

# Analyze performance by lead time
observations = continuous_observations()

# Calculate MSE for each lead time
mse_by_lead = mse(lead_forecast, observations,
                  reduce_dims=["time", "station"])

print(f"MSE by lead time:")
for lead in mse_by_lead.lead_time.values:
    lead_mse = mse_by_lead.sel(lead_time=lead)
    print(f"  Lead {lead:2d}: MSE = {lead_mse.values:.3f}")

from scores.sample_data import cdf_forecast, cdf_observations
from scores.probability import crps_cdf

# Generate CDF forecast data
cdf_fcst = cdf_forecast()
cdf_obs = cdf_observations()

print(f"CDF forecast dimensions: {cdf_fcst.dims}")
print(f"CDF thresholds: {len(cdf_fcst.threshold)} points")
print(f"Threshold range: [{cdf_fcst.threshold.min().values:.1f}, {cdf_fcst.threshold.max().values:.1f}]")

# Verify CDF properties
print(f"CDF starts at: {cdf_fcst.isel(threshold=0).values.mean():.3f}")
print(f"CDF ends at: {cdf_fcst.isel(threshold=-1).values.mean():.3f}")

# Calculate CRPS for CDF forecasts
crps_score = crps_cdf(cdf_fcst, cdf_obs, threshold_dim="threshold")
print(f"CRPS for CDF forecast: {crps_score.values:.3f}")

import matplotlib.pyplot as plt
import numpy as np

# Generate various sample data
simple_fcst = simple_forecast()
simple_obs = simple_observations()
cont_fcst = continuous_forecast()
cont_obs = continuous_observations()

# Create comparison plots
fig, ((ax1, ax2), (ax3, ax4)) = plt.subplots(2, 2, figsize=(12, 10))

# Simple data scatter plot
ax1.scatter(simple_fcst, simple_obs, alpha=0.6)
ax1.plot([simple_fcst.min(), simple_fcst.max()],
         [simple_fcst.min(), simple_fcst.max()], 'r--', alpha=0.7)
ax1.set_xlabel('Simple Forecast')
ax1.set_ylabel('Simple Observation')
ax1.set_title('Simple Data Scatter Plot')
ax1.grid(True)

# Time series plot (if time dimension exists)
if 'time' in cont_fcst.dims:
    time_slice = cont_fcst.isel(station=0) if 'station' in cont_fcst.dims else cont_fcst
    obs_slice = cont_obs.isel(station=0) if 'station' in cont_obs.dims else cont_obs

    ax2.plot(time_slice, label='Forecast', alpha=0.7)
    ax2.plot(obs_slice, label='Observation', alpha=0.7)
    ax2.set_xlabel('Time Index')
    ax2.set_ylabel('Value')
    ax2.set_title('Time Series Comparison')
    ax2.legend()
    ax2.grid(True)

# Distribution comparison
ax3.hist(simple_fcst, bins=20, alpha=0.6, label='Forecast', density=True)
ax3.hist(simple_obs, bins=20, alpha=0.6, label='Observation', density=True)
ax3.set_xlabel('Value')
ax3.set_ylabel('Density')
ax3.set_title('Distribution Comparison')
ax3.legend()
ax3.grid(True)

# Error analysis
errors = simple_fcst - simple_obs
ax4.hist(errors, bins=20, alpha=0.7, color='green')
ax4.axvline(0, color='red', linestyle='--', alpha=0.7)
ax4.set_xlabel('Forecast Error')
ax4.set_ylabel('Frequency')
ax4.set_title('Error Distribution')
ax4.grid(True)

plt.tight_layout()
plt.show()

# Summary statistics
print(f"\nData Summary Statistics:")
print(f"Simple forecast: mean={simple_fcst.mean().values:.3f}, std={simple_fcst.std().values:.3f}")
print(f"Simple observation: mean={simple_obs.mean().values:.3f}, std={simple_obs.std().values:.3f}")
print(f"Error statistics: mean={errors.mean().values:.3f}, std={errors.std().values:.3f}")

# Use sample data as templates for custom data generation
template_fcst = continuous_forecast()
template_obs = continuous_observations()

# Create custom data with similar structure but different values
custom_forecast = template_fcst.copy()
custom_forecast.values = np.random.normal(
    template_fcst.mean(),
    template_fcst.std() * 1.2,  # 20% more variable
    template_fcst.shape
)

custom_observations = template_obs.copy()
custom_observations.values = np.random.normal(
    template_obs.mean(),
    template_obs.std(),
    template_obs.shape
)

# Verify custom data maintains structure
print(f"Custom data verification:")
print(f"Dimensions match: {custom_forecast.dims == template_fcst.dims}")
print(f"Coordinates match: {list(custom_forecast.coords.keys()) == list(template_fcst.coords.keys())}")

# Score custom data
custom_mse = mse(custom_forecast, custom_observations)
print(f"Custom data MSE: {custom_mse.values:.3f}")

# Generate multiple datasets for ensemble analysis
n_datasets = 10
forecast_ensemble = []
observation_ensemble = []

for i in range(n_datasets):
    # Add some random seed variation
    np.random.seed(42 + i)
    fcst = continuous_forecast()
    obs = continuous_observations()

    forecast_ensemble.append(fcst)
    observation_ensemble.append(obs)

# Analyze ensemble statistics
ensemble_mses = [mse(f, o).values for f, o in zip(forecast_ensemble, observation_ensemble)]

print(f"Ensemble MSE Statistics:")
print(f"Mean MSE: {np.mean(ensemble_mses):.3f}")
print(f"MSE Range: [{np.min(ensemble_mses):.3f}, {np.max(ensemble_mses):.3f}]")
print(f"MSE Std Dev: {np.std(ensemble_mses):.3f}")