メニュー
Use this Skill for climate trend analysis: Mann-Kendall test, Sen's slope, extreme indices (RX1day, R10mm), percentile thresholds, IPCC-style figures.
デフォルトでは、最初にソースを確認する Prompt が選択されています。直接コマンドに切り替えるか、ローカルコピーをダウンロードすることもできます。
インストールを決める前に、SKILL.md と SkillsMP に表示されている付属ファイルをお読みください。
SOC 職業分類に基づく
Codex または Claude でインストール この Prompt をコピーして Codex、Claude、または他のアシスタントに貼り付けると、Skill ページを確認してインストールできます。
直接コマンドでは確認用 Prompt が省略されます。実行前にソースを確認してください。
npx skills add https://github.com/xjtulyc/awesome-rosetta-skills --skill climate-trendsコマンドは1行のまま表示されます。コピー前に横へスクロールして全体を確認してください。
ローカルで確認しますか？SkillsMP が現在取得できるファイルをダウンロードできます。
SKILL.md を表示中
Astronomical data analysis with astropy and astroquery — FITS I/O, WCS transforms, catalog cross-matching, aperture photometry, and CMB power spectra.
Seismological data analysis with ObsPy — FDSN waveform download, response removal, phase picking, moment tensor inversion, and seismicity mapping.
Download and analyze oceanographic data from Copernicus Marine Service and Argo floats using copernicusmarine, gsw, and xarray.
| name | climate-trends |
| description | Use this Skill for climate trend analysis: Mann-Kendall test, Sen's slope, extreme indices (RX1day, R10mm), percentile thresholds, IPCC-style figures. |
| tags | ["earth-science","climate","trend-analysis","xclim","extreme-events"] |
| version | 1.0.0 |
| authors | [{"name":"Rosetta Skills Contributors","github":"@xjtulyc"}] |
| license | MIT |
| platforms | ["claude-code","codex","gemini-cli","cursor"] |
| dependencies | {"python":["xclim>=0.45","scipy>=1.11","matplotlib>=3.7","numpy>=1.24","xarray>=2023.6","pandas>=2.0","cartopy>=0.22"]} |
| last_updated | 2026-03-17 |
| status | stable |
One-line summary: Compute climate trends (Mann-Kendall, Sen's slope), extreme precipitation/temperature indices, and produce IPCC AR6-style figures using xclim, scipy, and xarray.
Trigger keywords: Mann-Kendall trend test, Sen's slope, climate extremes, ETCCDI, RX1day, R10mm, TX90p, precipitation trend, temperature trend, IPCC figure
Non-parametric test for monotonic trends in time series:
$$ S = \sum_{i=1}^{n-1} \sum_{j=i+1}^{n} \text{sgn}(x_j - x_i) $$
Under $H_0$ (no trend), $S \sim N(0, V(S))$. Reject $H_0$ if $|Z| > z_{\alpha/2}$.
Robust, non-parametric trend magnitude:
$$ \hat{\beta} = \text{median}\left{ \frac{x_j - x_i}{j - i} : i < j \right} $$
| Index | Definition | Unit |
|---|---|---|
| RX1day | Max 1-day precipitation | mm |
| RX5day | Max 5-day precipitation | mm |
| R10mm | Annual days with precip ≥ 10mm | days |
| R20mm | Annual days with precip ≥ 20mm | days |
| TX90p | Days with Tmax > 90th percentile | % |
| TN10p | Days with Tmin < 10th percentile | % |
| WSDI | Warm spell duration index | days |
| CDD |
| Consecutive dry days |
| days |
pip install xclim>=0.45 scipy>=1.11 matplotlib>=3.7 numpy>=1.24 \
xarray>=2023.6 pandas>=2.0 cartopy>=0.22
import xclim
import xarray as xr
import scipy
import numpy as np
print(f"xclim: {xclim.__version__}")
print(f"xarray: {xr.__version__}")
print(f"scipy: {scipy.__version__}")
# Expected: xclim 0.45+, xarray 2023.6+
import numpy as np
import xarray as xr
import pandas as pd
import matplotlib.pyplot as plt
def create_synthetic_climate_dataset(n_years=60, seed=42):
"""
Create a synthetic daily climate dataset (temperature + precipitation).
Replace with: xr.open_dataset('era5_data.nc') or similar.
"""
rng = np.random.default_rng(seed)
dates = pd.date_range("1965-01-01", periods=n_years*365, freq="D")
years = dates.year
doy = dates.day_of_year
# Temperature with warming trend + seasonality
trend_temp = 0.02 * (years - 1965) # 0.02°C/year
seasonal_temp = 15 * np.sin(2 * np.pi * (doy - 90) / 365)
noise_temp = rng.normal(0, 2, len(dates))
tas = 15 + trend_temp + seasonal_temp + noise_temp # °C
# Precipitation with slight increasing trend + seasonality
trend_pr = 0.5 * (years - 1965) / n_years # 0.5mm/day trend over all years
seasonal_pr = 2 + 2 * np.cos(2 * np.pi * (doy - 200) / 365)
pr = np.maximum(0, seasonal_pr + trend_pr + rng.exponential(3, len(dates)) - 3)
ds = xr.Dataset({
"tas": xr.DataArray(tas, dims=["time"],
attrs={"units": "degC", "long_name": "Near-surface air temperature"}),
"pr": xr.DataArray(pr, dims=["time"],
attrs={"units": "mm/day", "long_name": "Precipitation"}),
}, coords={"time": dates})
return ds
ds = create_synthetic_climate_dataset()
print(f"Dataset: {len(ds.time)} days ({ds.time.dt.year.min().item()}–{ds.time.dt.year.max().item()})")
print(f"Temperature range: {float(ds.tas.min()):.1f} – {float(ds.tas.max()):.1f} °C")
print(f"Max daily precip: {float(ds.pr.max()):.1f} mm/day")
# Annual aggregates
tas_annual = ds.tas.resample(time="YE").mean()
pr_annual = ds.pr.resample(time="YE").sum() # annual total
years = tas_annual.time.dt.year.values
print(f"\nAnnual temp: {float(tas_annual.mean()):.2f} ± {float(tas_annual.std()):.2f} °C")
print(f"Annual precip: {float(pr_annual.mean()):.0f} ± {float(pr_annual.std()):.0f} mm/year")
import numpy as np
import scipy.stats as stats
from scipy.stats import kendalltau
import matplotlib.pyplot as plt
def mann_kendall_test(series):
"""
Mann-Kendall trend test and Sen's slope estimator.
Returns
-------
dict with: tau, p_value, sen_slope, confidence_interval
"""
x = np.asarray(series)
n = len(x)
# Kendall's tau
tau, p_value = kendalltau(np.arange(n), x)
# Sen's slope
slopes = []
for i in range(n):
for j in range(i+1, n):
slopes.append((x[j] - x[i]) / (j - i))
sen_slope = np.median(slopes)
# 95% CI for Sen's slope
C_alpha = 1.96 * np.sqrt(n*(n-1)*(2*n+5)/18)
M1 = int((len(slopes) - C_alpha) / 2)
M2 = int((len(slopes) + C_alpha) / 2)
slopes_sorted = sorted(slopes)
ci_lo = slopes_sorted[max(0, M1)]
ci_hi = slopes_sorted[min(len(slopes)-1, M2)]
return {
"tau": tau,
"p_value": p_value,
"sen_slope": sen_slope,
"ci_95": (ci_lo, ci_hi),
"significant": p_value < 0.05,
}
# Compute trends
tas_annual_np = float(tas_annual.mean()) + np.arange(len(tas_annual)) * 0.025 + \
np.random.normal(0, 0.5, len(tas_annual)) # example annual values
result_temp = mann_kendall_test(tas_annual_np)
print("Temperature trend analysis:")
print(f" Kendall τ = {result_temp['tau']:.4f}, p = {result_temp['p_value']:.4f}")
print(f" Sen's slope = {result_temp['sen_slope']:.4f} °C/year (95% CI: {result_temp['ci_95'][0]:.4f} – {result_temp['ci_95'][1]:.4f})")
print(f" Significant trend: {result_temp['significant']}")
years_arr = np.arange(1965, 1965 + len(tas_annual_np))
trend_line = result_temp['sen_slope'] * (years_arr - years_arr.mean()) + np.mean(tas_annual_np)
fig, ax = plt.subplots(figsize=(10, 5))
ax.bar(years_arr, tas_annual_np, color="salmon", alpha=0.6, label="Annual mean T")
ax.plot(years_arr, trend_line, 'r-', linewidth=2.5,
label=f"Sen's slope: {result_temp['sen_slope']*10:.2f} °C/decade")
ax.set_xlabel("Year"); ax.set_ylabel("Temperature (°C)")
ax.set_title(f"Annual Temperature Trend (τ={result_temp['tau']:.3f}, p={result_temp['p_value']:.3f})")
ax.legend()
plt.tight_layout()
plt.savefig("temperature_trend.png", dpi=150)
plt.show()
import xclim
import xarray as xr
import pandas as pd
import matplotlib.pyplot as plt
import numpy as np
# Create dataset with xclim-compatible attributes
ds = create_synthetic_climate_dataset()
# Convert units for xclim (expects Kelvin for temperature, kg/m2/s for precip)
tas_K = ds.tas + 273.15
tas_K.attrs["units"] = "K"
pr_kgm2s = ds.pr / 86400 # mm/day → kg/m2/s
pr_kgm2s.attrs["units"] = "kg m-2 s-1"
print("Computing ETCCDI indices...")
# Maximum 1-day precipitation (RX1day)
rx1day = xclim.indices.max_1day_precipitation_amount(pr=pr_kgm2s, freq="YS")
print(f"RX1day range: {float(rx1day.min()):.1f} – {float(rx1day.max()):.1f} mm")
# Days with precipitation >= 10mm
r10mm = xclim.indices.days_over_precip_thresh(
pr=pr_kgm2s,
per=xr.DataArray([10/86400], dims=["quantile"],
attrs={"units": "kg m-2 s-1"}),
freq="YS"
)
print(f"R10mm (days ≥ 10mm/day): {float(r10mm.mean()):.1f} days/year")
# Warm days (Tmax > 90th percentile)
# Compute 90th percentile climatology first
tasmax = tas_K + np.abs(np.random.normal(3, 1, len(tas_K))) # simulate Tmax
tasmax_da = xr.DataArray(tasmax, dims=["time"],
coords={"time": ds.time},
attrs={"units": "K"})
t90p = xclim.indices.tg90p(tas=tas_K, t90=tas_K.quantile(0.90).values, freq="YS")
fig, axes = plt.subplots(2, 2, figsize=(12, 8))
time_years = rx1day.time.dt.year.values
for ax, data, label in zip(axes.flat,
[rx1day*86400, r10mm, t90p, (tas_K - 273.15).resample(time="YE").mean()],
["RX1day (mm)", "R10mm (days)", "TX90p (%)", "Annual mean T (°C)"]):
arr = data.values.flatten()
ax.bar(time_years[:len(arr)], arr, alpha=0.7, color="steelblue")
# Add trend line
if len(arr) > 3:
slope, intercept, *_ = np.polyfit(range(len(arr)), arr, 1), None, None
slope = np.polyfit(range(len(arr)), arr, 1)[0]
trend = np.polyfit(range(len(arr)), arr, 1)
ax.plot(time_years[:len(arr)], np.polyval(trend, range(len(arr))),
'r--', linewidth=1.5, label=f"Trend: {slope:.3f}/yr")
ax.set_ylabel(label)
ax.legend(fontsize=8)
plt.suptitle("ETCCDI Climate Extreme Indices (1965–2024)")
plt.tight_layout()
plt.savefig("etccdi_indices.png", dpi=150)
plt.show()
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.colors as mcolors
def warming_stripes(annual_temps, start_year=1965, reference_period=(1981, 2010)):
"""
Create IPCC AR6-style warming stripes.
"""
temps = np.array(annual_temps)
years = np.arange(start_year, start_year + len(temps))
# Compute anomaly relative to reference period
ref_mask = (years >= reference_period[0]) & (years <= reference_period[1])
ref_mean = temps[ref_mask].mean()
anomaly = temps - ref_mean
# Colormap: blue (cool) → white (neutral) → red (warm)
cmap = mcolors.LinearSegmentedColormap.from_list(
"warming_stripes",
["#08306b", "#2171b5", "#6baed6", "#c6dbef", "#ffffff",
"#fcbba1", "#fb6a4a", "#cb181d", "#67000d"]
)
vmin, vmax = -np.abs(anomaly).max(), np.abs(anomaly).max()
norm = mcolors.TwoSlopeNorm(vmin=vmin, vcenter=0, vmax=vmax)
fig, ax = plt.subplots(figsize=(12, 3))
for i, (year, anom) in enumerate(zip(years, anomaly)):
ax.axvspan(year - 0.5, year + 0.5, color=cmap(norm(anom)))
ax.set_xlim(years[0] - 0.5, years[-1] + 0.5)
ax.axis("off")
ax.set_title(f"Temperature Anomaly ({start_year}–{start_year+len(temps)-1})",
fontsize=14, fontweight="bold", pad=10)
sm = plt.cm.ScalarMappable(cmap=cmap, norm=norm)
cbar = plt.colorbar(sm, ax=ax, orientation='horizontal',
fraction=0.04, pad=0.01, shrink=0.5)
cbar.set_label("Temperature anomaly (°C)", fontsize=10)
plt.tight_layout()
plt.savefig("warming_stripes.png", dpi=200, bbox_inches="tight")
plt.show()
return anomaly
rng = np.random.default_rng(42)
n = 60
simulated_temps = (15 + 0.025 * np.arange(n) +
rng.normal(0, 0.5, n))
anomalies = warming_stripes(simulated_temps)
print(f"Max anomaly: {anomalies.max():.2f} °C")
xclim.core.utils.ValidationError: unitsCause: Input DataArray units attribute missing or wrong.
Fix:
# Always set units attribute before calling xclim
pr_kgm2s.attrs["units"] = "kg m-2 s-1"
tas_K.attrs["units"] = "K"
Cause: Short time series (< 10 points) or too much variability.
Fix:
# Use scipy's kendalltau with method='auto' for small samples
tau, p = stats.kendalltau(x, y, method='auto')
| Package | Tested versions | Known issues |
|---|---|---|
| xclim | 0.45, 0.48 | Index API may change between major versions |
| xarray | 2023.6, 2024.3 | None |
# =============================================
# Extreme precipitation frequency analysis
# =============================================
import numpy as np, pandas as pd, matplotlib.pyplot as plt
from scipy import stats
rng = np.random.default_rng(42)
n_years = 60
years = np.arange(1965, 1965 + n_years)
# Simulate annual maximum 1-day precipitation (GEV distribution)
# In practice: use RX1day from real data
annual_max = stats.genextreme.rvs(c=0.1, loc=50, scale=15, size=n_years, random_state=42)
annual_max += 0.3 * (years - 1965) / n_years * annual_max # slight trend
# GEV fit
shape, loc, scale = stats.genextreme.fit(annual_max)
print(f"GEV fit: shape={shape:.3f}, loc={loc:.1f}, scale={scale:.1f}")
# Return levels
return_periods = [2, 5, 10, 25, 50, 100]
exceedance_probs = [1/rp for rp in return_periods]
return_levels = stats.genextreme.isf(exceedance_probs, shape, loc, scale)
fig, axes = plt.subplots(1, 2, figsize=(12, 5))
axes[0].bar(years, annual_max, alpha=0.7, color="steelblue")
slope, intercept, *_ = stats.linregress(years, annual_max)
axes[0].plot(years, slope*years + intercept, 'r--', linewidth=2,
label=f"Trend: {slope:.2f} mm/yr")
axes[0].set_xlabel("Year"); axes[0].set_ylabel("Annual max 1-day precip (mm)")
axes[0].set_title("Annual Maximum Precipitation"); axes[0].legend()
axes[1].semilogx(return_periods, return_levels, 'bs-', linewidth=2)
axes[1].set_xlabel("Return period (years)")
axes[1].set_ylabel("Return level (mm)")
axes[1].set_title("Precipitation Frequency Curve (GEV)")
axes[1].grid(True, which="both", alpha=0.3)
plt.tight_layout()
plt.savefig("precip_frequency.png", dpi=150)
plt.show()
for rp, rl in zip(return_periods, return_levels):
print(f" {rp:3d}-year return level: {rl:.1f} mm")
Interpreting these results: A 100-year return level of X mm means that, on average, this precipitation amount is exceeded once every 100 years. Positive trend in annual maxima suggests increasing extreme precipitation.
Last updated: 2026-03-17 | Maintainer: @xjtulyc Issues: GitHub Issues