Environmental Compliance Analytics¶

# Databricks notebook source
# MAGIC %md
# MAGIC # Bronze → Silver: AQI Calculation from AQS Hourly Data
# MAGIC Calculates Air Quality Index from raw pollutant concentrations
# MAGIC using EPA breakpoint tables per 40 CFR Part 58 Appendix G.

# COMMAND ----------

from pyspark.sql import functions as F
from pyspark.sql.types import DoubleType, StringType, StructType, StructField

# EPA AQI Breakpoint Table (PM2.5 24-hour, µg/m³)
# Source: https://aqs.epa.gov/aqsweb/documents/codetables/aqi_breakpoints.html
aqi_breakpoints_pm25 = [
    (0.0, 12.0, 0, 50, "Good"),
    (12.1, 35.4, 51, 100, "Moderate"),
    (35.5, 55.4, 101, 150, "Unhealthy for Sensitive Groups"),
    (55.5, 150.4, 151, 200, "Unhealthy"),
    (150.5, 250.4, 201, 300, "Very Unhealthy"),
    (250.5, 500.4, 301, 500, "Hazardous"),
]

breakpoint_schema = StructType([
    StructField("conc_lo", DoubleType()),
    StructField("conc_hi", DoubleType()),
    StructField("aqi_lo", DoubleType()),
    StructField("aqi_hi", DoubleType()),
    StructField("category", StringType()),
])

df_bp = spark.createDataFrame(aqi_breakpoints_pm25, schema=breakpoint_schema)

# COMMAND ----------

# Read Bronze AQS hourly data
df_aqs = spark.read.format("delta").load(
    "Tables/bronze_epa_aqs_hourly"
)

# Calculate AQI using EPA linear interpolation formula:
# AQI = ((AQI_hi - AQI_lo) / (Conc_hi - Conc_lo)) * (Conc - Conc_lo) + AQI_lo
df_aqi = (
    df_aqs
    .filter(F.col("parameter_code") == "88101")  # PM2.5
    .crossJoin(df_bp)
    .filter(
        (F.col("sample_measurement") >= F.col("conc_lo")) &
        (F.col("sample_measurement") <= F.col("conc_hi"))
    )
    .withColumn(
        "aqi_value",
        ((F.col("aqi_hi") - F.col("aqi_lo")) /
         (F.col("conc_hi") - F.col("conc_lo"))) *
        (F.col("sample_measurement") - F.col("conc_lo")) +
        F.col("aqi_lo")
    )
    .withColumn("aqi_value", F.round(F.col("aqi_value")).cast("int"))
    .select(
        "state_code", "county_code", "site_number",
        "date_local", "time_local", "sample_measurement",
        "aqi_value", "category", "latitude", "longitude"
    )
)

# Write to Silver
df_aqi.write.format("delta").mode("overwrite").option(
    "overwriteSchema", "true"
).saveAsTable("lh_silver.silver_epa_aqi_calculated")

print(f"AQI records calculated: {df_aqi.count()}")

# Databricks notebook source
# MAGIC %md
# MAGIC # Silver → Gold: TRI Facility Release Trend Analysis
# MAGIC Identifies facilities with increasing toxic release patterns
# MAGIC over rolling 5-year windows for regulatory priority targeting.

# COMMAND ----------

from pyspark.sql import functions as F
from pyspark.sql.window import Window

# Read Silver TRI data (cleansed, deduplicated)
df_tri = spark.read.format("delta").load("Tables/silver_epa_tri_releases")

# Calculate facility-level annual release totals
df_facility_annual = (
    df_tri
    .groupBy("trifid", "facility_name", "state", "county",
             "naics_code", "reporting_year", "latitude", "longitude")
    .agg(
        F.sum("total_releases_lbs").alias("total_releases_lbs"),
        F.sum("fugitive_air_lbs").alias("fugitive_air_lbs"),
        F.sum("stack_air_lbs").alias("stack_air_lbs"),
        F.sum("water_discharge_lbs").alias("water_discharge_lbs"),
        F.sum("land_disposal_lbs").alias("land_disposal_lbs"),
        F.countDistinct("chemical_name").alias("distinct_chemicals"),
        F.collect_set("chemical_name").alias("chemicals_reported"),
    )
)

# Calculate 5-year rolling trend using linear regression slope
year_window = Window.partitionBy("trifid").orderBy("reporting_year")

df_trend = (
    df_facility_annual
    .withColumn("year_rank", F.row_number().over(year_window))
    .withColumn("releases_prev_year", F.lag("total_releases_lbs").over(year_window))
    .withColumn(
        "yoy_change_pct",
        F.when(
            F.col("releases_prev_year") > 0,
            ((F.col("total_releases_lbs") - F.col("releases_prev_year"))
             / F.col("releases_prev_year") * 100)
        ).otherwise(None)
    )
)

# Identify facilities with consistently increasing releases (3+ consecutive years)
increasing_window = Window.partitionBy("trifid").orderBy("reporting_year").rowsBetween(-2, 0)

df_flagged = (
    df_trend
    .withColumn("increasing_flag", F.when(F.col("yoy_change_pct") > 0, 1).otherwise(0))
    .withColumn("consecutive_increases", F.sum("increasing_flag").over(increasing_window))
    .withColumn(
        "risk_tier",
        F.when(F.col("consecutive_increases") >= 3, "HIGH")
        .when(F.col("consecutive_increases") == 2, "MEDIUM")
        .otherwise("LOW")
    )
)

# Write Gold aggregation
df_flagged.write.format("delta").mode("overwrite").option(
    "overwriteSchema", "true"
).saveAsTable("lh_gold.gold_epa_tri_facility_trends")

high_risk = df_flagged.filter(F.col("risk_tier") == "HIGH").select("trifid").distinct().count()
print(f"High-risk facilities (3+ consecutive years increasing): {high_risk}")

# Databricks notebook source
# MAGIC %md
# MAGIC # Gold: Multi-Program Compliance Risk Scoring
# MAGIC Calculates a composite 0-100 compliance risk score per facility
# MAGIC across CAA, CWA, RCRA, and SDWA enforcement programs.

# COMMAND ----------

from pyspark.sql import functions as F

# Read Silver ECHO compliance data
df_echo = spark.read.format("delta").load("Tables/silver_epa_echo_compliance")

# Define severity weights by violation type
severity_weights = {
    "HPV": 10.0,   # High Priority Violation (CAA)
    "SNC": 10.0,   # Significant Non-Compliance (CWA)
    "SNY": 8.0,    # Significant Non-Complier (RCRA)
    "SV":  6.0,    # Serious Violator (SDWA)
    "V":   3.0,    # Violation
    "N":   0.0,    # No Violation
}

# Calculate component scores
df_scored = (
    df_echo
    .withColumn(
        "severity_score",
        F.when(F.col("compliance_status") == "HPV", 10.0)
        .when(F.col("compliance_status") == "SNC", 10.0)
        .when(F.col("compliance_status") == "SNY", 8.0)
        .when(F.col("compliance_status") == "SV", 6.0)
        .when(F.col("compliance_status") == "V", 3.0)
        .otherwise(0.0)
    )
    .withColumn(
        "recency_weight",
        F.when(F.datediff(F.current_date(), F.col("violation_date")) <= 365, 1.0)
        .when(F.datediff(F.current_date(), F.col("violation_date")) <= 730, 0.7)
        .when(F.datediff(F.current_date(), F.col("violation_date")) <= 1095, 0.4)
        .otherwise(0.2)
    )
    .withColumn("weighted_severity", F.col("severity_score") * F.col("recency_weight"))
)

# Aggregate to facility level
df_facility_risk = (
    df_scored
    .groupBy("registry_id", "facility_name", "state", "city",
             "zip_code", "latitude", "longitude")
    .agg(
        F.sum("weighted_severity").alias("raw_risk_score"),
        F.count("*").alias("total_violations"),
        F.countDistinct("program_code").alias("programs_violated"),
        F.max("violation_date").alias("most_recent_violation"),
        F.sum(F.when(F.col("compliance_status").isin("HPV", "SNC"), 1).otherwise(0))
            .alias("high_priority_count"),
    )
    .withColumn(
        "composite_risk_score",
        F.least(
            F.lit(100),
            F.round(
                F.col("raw_risk_score") +
                (F.col("programs_violated") * 5) +
                (F.col("high_priority_count") * 3)
            )
        ).cast("int")
    )
    .withColumn(
        "risk_category",
        F.when(F.col("composite_risk_score") >= 75, "CRITICAL")
        .when(F.col("composite_risk_score") >= 50, "HIGH")
        .when(F.col("composite_risk_score") >= 25, "MODERATE")
        .otherwise("LOW")
    )
)

df_facility_risk.write.format("delta").mode("overwrite").option(
    "overwriteSchema", "true"
).saveAsTable("lh_gold.gold_epa_facility_risk_scores")

critical = df_facility_risk.filter(F.col("risk_category") == "CRITICAL").count()
print(f"Critical-risk facilities: {critical}")

# Cross-domain: USDA crop intensity vs. EPA water violations by county
df_cross_ag = (
    df_usda_crop_production
    .join(df_epa_water_violations, on=["state_fips", "county_fips"])
    .groupBy("state", "county", "reporting_year")
    .agg(
        F.sum("total_fertilizer_tons").alias("fertilizer_intensity"),
        F.sum("livestock_head_count").alias("livestock_density"),
        F.count("violation_id").alias("water_violations"),
        F.avg("nitrate_concentration_mg_l").alias("avg_nitrate_level"),
    )
    .withColumn(
        "runoff_risk_score",
        (F.col("fertilizer_intensity") * 0.4) +
        (F.col("livestock_density") * 0.3) +
        (F.col("water_violations") * 0.3)
    )
)

# Cross-domain: EPA contaminated facilities on DOI-managed lands
from pyspark.sql import functions as F

df_cross_lands = (
    df_epa_tri_facilities
    .join(
        df_doi_land_boundaries,
        (F.col("tri.latitude").between(F.col("land.lat_min"), F.col("land.lat_max"))) &
        (F.col("tri.longitude").between(F.col("land.lon_min"), F.col("land.lon_max")))
    )
    .groupBy("land_unit_name", "managing_bureau", "state")
    .agg(
        F.count("trifid").alias("tri_facilities_on_land"),
        F.sum("total_releases_lbs").alias("total_releases_on_land"),
        F.countDistinct("chemical_name").alias("distinct_chemicals"),
    )
    .orderBy(F.col("total_releases_on_land").desc())
)