مقدمة في تحليل البيانات الزراعية بـ Python

import pandas as pd
import numpy as np

df = pd.read_csv("datasets/kfs_farm_hourly_2023_2024.csv", parse_dates=["time"])
df = df.set_index("time").sort_index()        # DatetimeIndex -> time superpowers

df.head()            # first 5 rows
df.info()            # types + missing values
df.describe().T      # descriptive statistics for every column
df.loc["2024-07-15"] # all hours of one day (slicing by date string!)

stats = df[["air_temp_c", "air_humidity_pct", "soil_moisture_0_7cm", "et0_mm"]].describe().round(2)
print(stats)
print("hours above 40 °C:", (df["air_temp_c"] > 40).sum())
print("10th / 90th percentile of topsoil VWC:", df["soil_moisture_0_7cm"].quantile([.1, .9]).round(3).tolist())

monthly = df.resample("MS").agg({"air_temp_c": "mean", "precipitation_mm": "sum", "et0_mm": "sum"}).round(1)
monthly["water_balance_mm"] = monthly["precipitation_mm"] - monthly["et0_mm"]   # negative = irrigation deficit
print(monthly.loc["2023"])
print("Annual ET0 2023 (mm):", round(df.loc["2023", "et0_mm"].sum()))
print("Annual rain 2023 (mm):", round(df.loc["2023", "precipitation_mm"].sum(), 1))

summer = df.loc["2024-06":"2024-08"]
diurnal = summer.groupby(summer.index.hour)[["air_temp_c", "air_humidity_pct", "et0_mm"]].mean().round(2)
print(diurnal)
print("hottest hour:", diurnal["air_temp_c"].idxmax(), "| peak ET0 hour:", diurnal["et0_mm"].idxmax())
# agronomic reading: irrigate early morning (lowest ET0 + lowest temperature) to cut evaporation losses

dry = df["soil_moisture_0_7cm"] < 0.110            # boolean mask (threshold = 20th percentile of this cell)
block = (dry != dry.shift()).cumsum()              # new id every time the state flips
periods = df[dry].groupby(block[dry]).agg(
    start=("soil_moisture_0_7cm", lambda s: s.index.min()),
    hours=("soil_moisture_0_7cm", "size"))
print(periods[periods.hours >= 48].sort_values("hours", ascending=False).head())

"""Session 13 - Agricultural data analysis with Pandas & NumPy.

Goal: describe the climate/soil of the KFS farm and detect drought periods.
Run:  python s13_explore_farm_data.py
Outputs (outputs/): s13_weekly_report.csv, s13_drought_periods.csv, s13_overview.png
"""
import matplotlib
matplotlib.use("Agg")
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd

pd.set_option("display.width", 160)
pd.set_option("display.max_columns", 20)

from common import OUT, RAINFED_DRY as REFILL_POINT, load_hourly

df = load_hourly()
print("Rows:", len(df), "| from", df.index.min(), "to", df.index.max())
print("\n== Descriptive statistics ==")
print(df.describe().T[["mean", "std", "min", "max"]].round(2))

# ---- 1. Grouping by time: monthly climate profile
monthly = df.resample("MS").agg({"air_temp_c": "mean", "air_humidity_pct": "mean",
                                 "precipitation_mm": "sum", "et0_mm": "sum",
                                 "soil_moisture_0_7cm": "mean"}).round(2)
monthly["water_balance_mm"] = (monthly["precipitation_mm"] - monthly["et0_mm"]).round(1)
print("\n== Monthly profile ==")
print(monthly.head(12))

# ---- 2. Diurnal cycle: when is it hottest? (helps choose irrigation time)
diurnal = df.groupby(df.index.hour)[["air_temp_c", "air_humidity_pct", "et0_mm"]].mean().round(2)
print("\n== Hour of hottest air (mean):", int(diurnal["air_temp_c"].idxmax()), ":00")

# ---- 3. Weekly report for one season (summer 2024)
summer = df.loc["2024-06-01":"2024-08-31"]
weekly = summer.resample("W").agg(temp_mean=("air_temp_c", "mean"), temp_max=("air_temp_c", "max"),
                                  rain=("precipitation_mm", "sum"), et0=("et0_mm", "sum"),
                                  soil_mean=("soil_moisture_0_7cm", "mean"),
                                  soil_min=("soil_moisture_0_7cm", "min")).round(3)
weekly.to_csv(OUT / "s13_weekly_report.csv")

# ---- 4. Drought periods: consecutive hours with soil moisture below the refill point
dry = df["soil_moisture_0_7cm"] < REFILL_POINT
block_id = (dry != dry.shift()).cumsum()
periods = (df[dry].groupby(block_id[dry])
           .agg(start=("soil_moisture_0_7cm", lambda s: s.index.min()),
                end=("soil_moisture_0_7cm", lambda s: s.index.max()),
                hours=("soil_moisture_0_7cm", "size"),
                min_vwc=("soil_moisture_0_7cm", "min")))
periods = periods[periods["hours"] >= 48].sort_values("hours", ascending=False)
periods.to_csv(OUT / "s13_drought_periods.csv", index=False)
print(f"\n== Drought periods (>=48 h below {REFILL_POINT} m3/m3): {len(periods)}")
print(periods.head(5).to_string(index=False))

# ---- 5. Figure
fig, ax = plt.subplots(3, 1, figsize=(11, 8), sharex=True)
daily = df.resample("D").mean(numeric_only=True)
ax[0].plot(daily.index, daily["air_temp_c"], color="#e07a1f"); ax[0].set_ylabel("Air temp °C")
ax[1].bar(daily.index, df["precipitation_mm"].resample("D").sum(), color="#2a6fdb"); ax[1].set_ylabel("Rain mm")
ax[2].plot(daily.index, daily["soil_moisture_0_7cm"], color="#13787b")
ax[2].axhline(REFILL_POINT, ls="--", color="crimson", label="refill point"); ax[2].legend()
ax[2].set_ylabel("Soil VWC m³/m³")
fig.suptitle("KFS farm - climate & soil overview (2023-2024)")
fig.tight_layout(); fig.savefig(OUT / "s13_overview.png", dpi=120)
print("\nSaved outputs to", OUT)
print("Correlation soil moisture vs ET0 (daily):",
      np.round(daily["soil_moisture_0_7cm"].corr(df["et0_mm"].resample("D").sum()), 3))

import sqlite3
con = sqlite3.connect("smart-farm/farm.db")
nodes = pd.read_sql("SELECT node_id, ts, soil_moisture_pct, air_temp_c FROM readings", con, parse_dates=["ts"])
hourly = nodes.set_index("ts").groupby("node_id")["soil_moisture_pct"].resample("h").mean().unstack(0)
report = pd.DataFrame({"mean": hourly.mean(), "min": hourly.min(), "hours_below_30": (hourly < 30).sum()}).round(1)
print(report)