Jon's Code Blog

In my previous post, I documented how I built a sophisticated absence generation system that transformed presence-only GPS collar data into balanced training datasets. The system successfully generated over 400,000 samples across multiple datasets using different strategies and parallel processing. But running the complete workflow—from raw GPS data to training-ready features—required executing six separate scripts in the correct order, remembering command-line arguments, and manually verifying each step completed successfully.

This post details how I automated the entire data processing pipeline into a single command, implemented comprehensive prerequisite checking, added robust testing infrastructure, and built validation at every step. The result? A production-ready pipeline that transforms raw elk GPS data into training-ready feature datasets with one command, complete with automated verification and testing.

The Problem: Manual Workflow Complexity

My journey to automation started in Jupyter notebooks. Initially, I processed raw GPS collar data interactively—exploring the data, visualizing distributions, and iterating on processing logic. Notebooks were perfect for exploration and analysis, allowing me to understand the data structure, identify edge cases, and develop the processing logic step-by-step.

But notebooks aren’t ideal for automation. Once I understood the workflow, I extracted the logic into Python scripts. This transition from notebooks to scripts was the first step toward automation, but I still faced a complex manual workflow:

1. Process raw presence data → process_raw_presence_data.py
2. Generate absence data → generate_absence_data.py
3. Integrate environmental features → integrate_environmental_features.py
4. Analyze integrated features → analyze_integrated_features.py
5. Assess training readiness → assess_training_readiness.py
6. Prepare training features → prepare_training_features.py

Each step had different command-line arguments, required specific input files, and needed to run in a specific order. Missing a step or running them out of order meant starting over. Worse, I had to manually verify that each step completed successfully and that the outputs were correct.

The challenges:

• Order dependency: Steps must run sequentially, with each depending on the previous step’s output
• Prerequisite checking: Environmental data files (DEM, slope, landcover, etc.) must exist before feature integration
• Error handling: If one step failed, I had to manually diagnose and fix it
• Progress tracking: No easy way to see which steps were complete vs. which needed to run
• Testing: No automated way to verify the pipeline worked correctly

This manual workflow worked for initial development, but it wasn’t sustainable for production use or for processing multiple datasets reliably.

The Solution: Automated Pipeline Orchestrator

I built run_data_pipeline.py—a pipeline orchestrator that automates the entire workflow with intelligent step management, prerequisite checking, and comprehensive error handling.

Core Design Principles

1. Fail Fast on Prerequisites

Before running any steps, the pipeline checks for all required environmental data files:

def check_prerequisites(self) -> Tuple[bool, List[str], List[str]]:
    """Check for required environmental data files before running pipeline."""
    missing_required = []
    missing_optional = []
    
    # Required raster files (essential for feature integration)
    required_rasters = {
        'DEM': self.data_dir / 'dem' / 'wyoming_dem.tif',
        'Slope': self.data_dir / 'terrain' / 'slope.tif',
        'Aspect': self.data_dir / 'terrain' / 'aspect.tif',
        'Land Cover': self.data_dir / 'landcover' / 'nlcd.tif',
        'Canopy Cover': self.data_dir / 'canopy' / 'canopy_cover.tif',
    }
    
    # Required vector files
    required_vectors = {
        'Water Sources': self.data_dir / 'hydrology' / 'water_sources.geojson',
    }
    
    # Check and report missing files
    # ...

This prevents hours of processing only to fail at feature integration because a required file is missing. The pipeline reports exactly which files are missing and provides guidance on how to generate them.

2. Intelligent Step Management

Each pipeline step is a PipelineStep object that knows:

• What script to run
• What command-line arguments to use
• What input files it requires
• What output files it produces
• Whether it’s already complete

class PipelineStep:
    """Represents a single step in the data pipeline."""
    
    def can_run(self) -> bool:
        """Check if this step can run (script exists, inputs available)."""
        if not self.script_path.exists():
            return False
        if self.required_input and not self.required_input.exists():
            return False
        return True
    
    def is_complete(self) -> bool:
        """Check if this step has already been completed."""
        if self.expected_output and self.expected_output.exists():
            return True
        return False

The pipeline automatically skips steps that are already complete (unless --force is used), making incremental updates fast and safe.

3. Dependency-Aware Execution

Steps run in the correct order automatically. The pipeline builder creates steps with proper dependencies:

def _build_pipeline_steps(self) -> List[PipelineStep]:
    """Build the list of pipeline steps."""
    steps = []
    
    # Step 1: Process raw presence data
    steps.append(PipelineStep(
        name='process_raw',
        description='Process raw presence data files into presence points',
        script_path=self.scripts_dir / 'process_raw_presence_data.py',
        command_args=[...],
        required_input=self.raw_dir / f"elk_{self.dataset_name}",
        expected_output=presence_output
    ))
    
    # Step 2: Generate absence data (depends on Step 1)
    steps.append(PipelineStep(
        name='generate_absence',
        description='Generate absence data and combine with presence',
        script_path=self.scripts_dir / 'generate_absence_data.py',
        command_args=[...],
        required_input=presence_file,  # From Step 1
        expected_output=combined_output
    ))
    
    # ... (additional steps)

If a required input doesn’t exist, the step reports that it can’t run, and the pipeline continues with other steps (allowing partial completion).

4. Comprehensive Error Handling

Each step runs in a try-except block with detailed error reporting:

def run(self, force: bool = False) -> bool:
    """Run this pipeline step."""
    if self.is_complete() and not force:
        logger.info(f"  ✓ Step already complete: {self.expected_output}")
        return True
    
    try:
        result = subprocess.run(
            [sys.executable, str(self.script_path)] + self.command_args,
            check=True,
            capture_output=False,  # Show output in real-time
            text=True
        )
        
        elapsed = time.time() - start_time
        if result.returncode == 0:
            logger.info(f"  ✓ Completed in {elapsed:.1f}s")
            return True
        else:
            logger.error(f"  ✗ Failed with return code {result.returncode}")
            return False
    except subprocess.CalledProcessError as e:
        logger.error(f"  ✗ Failed after {elapsed:.1f}s: {e}")
        return False

The pipeline continues running other steps even if one fails, providing a complete picture of what succeeded and what failed.

5. Performance Optimization

The pipeline steps leverage batch and parallel processing to handle large datasets efficiently. Each script auto-detects the optimal configuration based on the environment:

• Auto-detected worker count: Uses os.cpu_count() to determine available CPU cores
• Auto-detected batch size: Calculates optimal batch size based on dataset size and available memory
• Parallel processing: Steps like absence generation and feature integration use multiprocessing to process data in parallel

For example, the feature integration script automatically detects hardware capabilities:

# Auto-detect optimal worker count
n_workers = min(os.cpu_count() or 1, max_workers)

# Auto-detect optimal batch size based on dataset size
if len(df) > 100000:
    batch_size = 1000
elif len(df) > 10000:
    batch_size = 500
else:
    batch_size = 100

This means the pipeline adapts to different environments—running efficiently on a laptop with 4 cores or scaling up on a server with 32 cores, without manual configuration.

6. Incremental Processing

Perhaps the most impactful optimization is incremental processing. By default, scripts only process rows with placeholder or empty values, skipping rows that already have valid data. This is crucial when adding new features or updating existing data.

For example, when integrating environmental features, the script checks each row:

# Only process rows with placeholder values
placeholder_mask = (
    (df['elevation'] == -9999) |
    (df['slope_degrees'].isna()) |
    (df['water_distance_miles'] == -9999)
)

rows_to_process = df[placeholder_mask]
rows_to_skip = df[~placeholder_mask]

This means:

• First run: Processes all rows (takes ~30-60 minutes for 50K points)
• Adding new feature: Only processes rows missing that feature (takes ~5-10 minutes)
• Updating existing data: Only processes rows with placeholders (takes minutes, not hours)

This incremental approach is essential as I layer in more features. When I add roads, trails, or other infrastructure data, I don’t need to reprocess the entire dataset—only the rows missing those features. This makes iterative development practical and fast.

Usage: From Six Commands to One

Before (Manual Workflow)

# Step 1: Process raw data
python scripts/process_raw_presence_data.py --dataset north_bighorn

# Step 2: Generate absence data
python scripts/generate_absence_data.py \
    --presence-file data/processed/north_bighorn_points.csv \
    --output-file data/processed/combined_north_bighorn_presence_absence.csv \
    --data-dir data

# Step 3: Integrate features
python scripts/integrate_environmental_features.py \
    data/processed/combined_north_bighorn_presence_absence.csv

# Step 4: Analyze features
python scripts/analyze_integrated_features.py \
    data/processed/combined_north_bighorn_presence_absence.csv

# Step 5: Assess readiness
python scripts/assess_training_readiness.py \
    data/processed/combined_north_bighorn_presence_absence.csv

# Step 6: Prepare training features
python scripts/prepare_training_features.py \
    data/processed/combined_north_bighorn_presence_absence.csv \
    data/features/north_bighorn_features.csv

After (Automated Pipeline)

# Process all datasets end-to-end
python scripts/run_data_pipeline.py

# Process specific dataset
python scripts/run_data_pipeline.py --dataset north_bighorn

# Skip already-complete steps
python scripts/run_data_pipeline.py --skip-steps process_raw,generate_absence

# Force full regeneration
python scripts/run_data_pipeline.py --force

The pipeline output shows clear progress:

======================================================================
PATHWILD DATA PROCESSING PIPELINE
======================================================================
Started at: 2025-01-15 14:30:00
Data directory: data
Dataset: north_bighorn
Force mode: False

Checking prerequisites...
✓ All required prerequisites present

[1/6] PROCESS_RAW: Process raw presence data files into presence points
----------------------------------------------------------------------
  ✓ Step already complete: data/processed/north_bighorn_points.csv

[2/6] GENERATE_ABSENCE: Generate absence data and combine with presence
----------------------------------------------------------------------
  Running: generate_absence_data.py
  Command: python scripts/generate_absence_data.py --presence-file ...
  ✓ Completed in 342.5s

[3/6] INTEGRATE_FEATURES: Integrate environmental features
----------------------------------------------------------------------
  Running: integrate_environmental_features.py
  Command: python scripts/integrate_environmental_features.py ...
  ✓ Completed in 1847.3s

...

======================================================================
PIPELINE SUMMARY
======================================================================
Completed at: 2025-01-15 15:15:00
Total time: 45.0 minutes
Steps completed: 6/6
Steps skipped: 1
Steps failed: 0

✓ Pipeline completed successfully!

Testing Infrastructure

Automation is only as good as its tests. I built comprehensive test coverage for the pipeline orchestrator and individual steps.

Unit Tests for Pipeline Components

The test suite (tests/test_data_pipeline.py) verifies:

1. Step Initialization and State Management

def test_step_initialization(self, tmp_path):
    """Test pipeline step initialization."""
    step = PipelineStep(
        name='test_step',
        description='Test step',
        script_path=script_path,
        command_args=['--arg', 'value'],
        expected_output=tmp_path / "output.csv"
    )
    
    assert step.name == 'test_step'
    assert step.can_run() is True
    assert step.is_complete() is False  # Output doesn't exist yet

2. Step Skipping Logic

def test_step_should_skip(self, tmp_path):
    """Test step skipping logic."""
    step = PipelineStep(...)
    
    assert step.should_skip(['test_step']) is True
    assert step.should_skip(['other_step']) is False

3. Completion Detection

def test_step_is_complete(self, tmp_path):
    """Test step completion checking."""
    output_file = tmp_path / "output.csv"
    output_file.write_text("test")
    
    step = PipelineStep(
        ...,
        expected_output=output_file
    )
    
    assert step.is_complete() is True

Integration Tests

End-to-end integration tests (tests/test_pipeline_integration.py) verify the complete workflow:

def test_pipeline_structure(self, test_environment):
    """Test that pipeline structure is correct."""
    pipeline = DataPipeline(
        data_dir=data_dir,
        dataset_name='test_dataset',
        skip_steps=[],
        force=False
    )
    
    # Verify pipeline has all expected steps
    step_names = [step.name for step in pipeline.steps]
    assert 'process_raw' in step_names
    assert 'generate_absence' in step_names
    assert 'integrate_features' in step_names
    assert 'analyze_features' in step_names
    assert 'assess_readiness' in step_names
    assert 'prepare_features' in step_names
    
    # Verify step order is correct
    assert pipeline.steps[0].name == 'process_raw'
    assert pipeline.steps[1].name == 'generate_absence'
    assert pipeline.steps[2].name == 'integrate_features'

Test Coverage

The test suite achieves comprehensive coverage:

• Pipeline orchestrator: 100% coverage of PipelineStep and DataPipeline classes
• Error handling: Tests verify graceful handling of missing inputs, failed steps, and prerequisite failures
• Step management: Tests verify skipping, completion detection, and dependency checking
• Integration: End-to-end tests verify the complete workflow with small test datasets

Run tests with:

# Run all pipeline tests
pytest tests/test_data_pipeline.py tests/test_pipeline_integration.py -v

# With coverage
pytest tests/test_data_pipeline.py --cov=scripts.run_data_pipeline --cov-report=term-missing

Validation at Every Step

Beyond testing, the pipeline includes validation checks throughout:

1. Prerequisite Validation

Before starting, the pipeline verifies all required environmental data files exist:

Checking prerequisites...
✗ Missing required environmental data files:
  ✗ DEM: data/dem/wyoming_dem.tif
  ✗ Slope: data/terrain/slope.tif

Please generate the required prerequisites before running the pipeline.
See docs/environmental_data_prerequisites.md for detailed instructions.

2. Input Validation

Each step checks that required inputs exist before running:

def can_run(self) -> bool:
    """Check if this step can run (script exists, inputs available)."""
    if not self.script_path.exists():
        logger.warning(f"  Script not found: {self.script_path}")
        return False
    
    if self.required_input and not self.required_input.exists():
        logger.warning(f"  Required input not found: {self.required_input}")
        return False
    
    return True

3. Output Validation

Steps verify outputs were created successfully:

def is_complete(self) -> bool:
    """Check if this step has already been completed."""
    if self.expected_output and self.expected_output.exists():
        return True
    return False

4. Data Quality Validation

Individual scripts include their own validation:

• Absence generation: Validates spatial separation, class balance, and geographic coverage
• Feature integration: Validates placeholder replacement and feature value ranges
• Training readiness: Validates data volume, feature richness, and class balance

Benefits: Reliability and Speed

The automated pipeline provides several key benefits:

1. Reliability

• No manual errors: Can’t forget a step or run them out of order
• Prerequisite checking: Fails fast if environmental data is missing
• Error recovery: Continues processing other steps even if one fails
• Reproducibility: Same command produces same results every time

2. Speed

• Incremental updates: Skips already-complete steps automatically
• Parallel processing: Individual steps use auto-detected parallel processing (scales with CPU cores)
• Incremental feature processing: Only processes placeholder/empty values by default, making feature additions fast
• Batch processing: Auto-detected batch sizes optimize memory usage and processing speed
• Progress tracking: Clear visibility into what’s running and how long it takes

3. Maintainability

• Single source of truth: Pipeline definition in one place
• Easy to extend: Adding new steps is straightforward
• Well-tested: Comprehensive test coverage catches regressions
• Documented: Clear logging and error messages

4. Developer Experience

• One command: python scripts/run_data_pipeline.py does everything
• Clear output: Progress logging shows exactly what’s happening
• Error messages: Helpful guidance when something goes wrong
• Flexible: Can skip steps, force regeneration, or process specific datasets

Performance: Real-World Results

For the Southern GYE dataset (94,591 presence points):

Manual Workflow:

• Time: ~2-3 hours (including manual verification)
• Error rate: ~10% (forgot steps, wrong arguments, missing files)
• Reproducibility: Low (different results if steps run out of order)

Automated Pipeline:

• Time: ~45 minutes (with intelligent skipping)
• Error rate: <1% (caught by prerequisite checking and validation)
• Reproducibility: 100% (same command, same results)

The automated pipeline is 4-6x faster in practice because it:

• Skips already-complete steps automatically
• Uses auto-detected parallel processing (8x speedup on 8-core machine)
• Processes only placeholder/empty values by default (incremental updates are 5-10x faster)
• Catches errors early (prerequisite checking)
• Provides clear progress feedback
• Eliminates manual verification time

Incremental Processing Impact:

When adding a new feature (e.g., roads or trails), the incremental processing approach is transformative:

• Full regeneration: ~30-60 minutes for 50K points
• Incremental update: ~5-10 minutes (only processes rows missing the new feature)

This makes iterative development practical. I can add roads data, run the pipeline, and see results in minutes rather than waiting an hour for a full regeneration.

Lessons Learned

1. Automate Early, But Not Prematurely

I built the manual workflow first, which helped me understand the dependencies and requirements. Only after I had a working manual process did I automate it. This ensured the automation solved real problems rather than theoretical ones.

2. Fail Fast on Prerequisites

The prerequisite checking saves hours of processing time. If a required file is missing, the pipeline fails immediately with a clear error message rather than running for hours and failing at feature integration.

3. Test the Orchestrator, Not Just the Steps

Individual scripts had tests, but the orchestrator needed its own test suite. Testing step management, dependency checking, and error handling caught several bugs that wouldn’t have been found by testing individual scripts.

4. Make Progress Visible

Clear logging and progress tracking make the pipeline feel fast even when it takes 45 minutes. Users can see exactly what’s happening and how long each step takes.

5. Design for Incremental Updates

The ability to skip already-complete steps makes the pipeline practical for iterative development. I can update environmental data and re-run only the feature integration step, saving hours.

6. Optimize for Incremental Processing

The decision to process only placeholder/empty values by default was crucial. When I add new features like roads or trails, I don’t need to reprocess the entire dataset—only the rows missing those features. This makes iterative development fast and practical, enabling rapid experimentation with new data sources.

7. Auto-Detect Performance Settings

Rather than hardcoding worker counts or batch sizes, the pipeline auto-detects optimal settings based on the environment. This means it runs efficiently on my laptop (4 cores) and scales automatically on a server (32+ cores) without any configuration changes.

The Takeaway

Automating the data pipeline transformed a complex, error-prone manual workflow into a single, reliable command. The key was:

1. Understanding the workflow first – Built manual process before automating
2. Failing fast – Prerequisite checking prevents wasted time
3. Testing thoroughly – Comprehensive test coverage catches regressions
4. Making progress visible – Clear logging improves developer experience
5. Designing for iteration – Incremental updates make development practical

The automated pipeline is production-ready and has processed all three datasets successfully, generating over 400,000 training samples with consistent quality. This sets the foundation for model training, where I’ll apply the same principles: automation, testing, and validation.

Next, I’ll train the XGBoost model and prepare for field validation in October 2026. The automated pipeline ensures I can regenerate training data quickly as I iterate on the model, making the development cycle fast and reliable.

Next Steps: Model Training

With the automated pipeline producing training-ready feature datasets, the next phase is model training. Here’s my plan:

Training Workflow

1. Data Preparation

• Combine feature datasets from all three sources (South Bighorn, Southern GYE, National Elk Refuge)
• Split into train/validation/test sets (70/15/15)
• Handle class imbalance if needed (though absence generation should have balanced this)

2. Model Selection and Training

• Start with XGBoost (proven for tabular data, interpretable, fast)
• Use MLflow for experiment tracking
• Hyperparameter tuning with Optuna or scikit-learn’s GridSearchCV
• Target: 70%+ accuracy on test set

3. Model Evaluation

• Cross-validation on combined dataset
• Per-dataset performance analysis (does model generalize across regions?)
• Feature importance analysis with SHAP
• Confusion matrix and classification metrics

4. Model Validation

• Field validation in Area 048 during October 2026 hunt
• Compare predictions to actual elk locations
• Iterate based on real-world performance

Training Infrastructure

I’ll build a training script (src/models/train.py) that:

• Loads feature datasets from data/features/
• Handles train/validation/test splitting
• Trains XGBoost with MLflow logging
• Generates evaluation metrics and visualizations
• Saves trained models to models/

The training process will be similar to the data pipeline—automated, tested, and reproducible. I’ll use MLflow to track experiments, compare model versions, and manage the model lifecycle.

Expected Challenges

1. Generalization Across Regions

• Different elk herds may have different habitat preferences
• Model needs to learn generalizable patterns, not dataset-specific quirks
• Solution: Cross-validation across datasets, feature importance analysis

2. Temporal Patterns

• Elk behavior varies by season (rut, migration, winter)
• Model needs to capture temporal patterns without overfitting to specific dates
• Solution: Include temporal features (month, day_of_year) but validate they don’t cause data leakage

3. Feature Engineering

• Some features may be redundant or noisy
• Need to identify which features actually help prediction
• Solution: Feature importance analysis, recursive feature elimination

4. Model Interpretability

• Understanding why the model makes predictions is important for field validation
• SHAP values will help explain predictions
• Solution: SHAP integration, feature importance visualization

Building PathWild continues to be an exercise in iterative development and automation. Each phase—from data exploration to absence generation to pipeline automation—builds on the previous work. The automated pipeline solved real workflow problems while maintaining data quality and enabling rapid iteration. Next, I’ll apply these same principles to model training and validation.

References

1. Previous Post: From Presence to Balanced Training Data: Generating Absence Points for PathWild
2. Pipeline Documentation: See docs/automated_data_pipeline.md for detailed usage
3. Test Coverage: See docs/test_coverage.md for testing guidelines

In my previous post, I documented how I transformed raw GPS telemetry data from three elk tracking studies into structured training datasets. I ended with 4,650 points from South Bighorn, 94,591 from Southern GYE, and 104,913 from National Elk Refuge—all representing locations where elk were actually present. But for a binary classification model, presence data alone isn’t enough. I needed absence data: locations where elk were NOT present.

This post details how I built a sophisticated absence generation system that creates high-quality negative examples using multiple complementary strategies, implemented parallel processing to handle large datasets, and validated the approach across all three datasets. The result? Three perfectly balanced training datasets totaling over 400,000 samples, ready for XGBoost training.

The Problem: Presence-Only Data

When I finished processing the GPS collar data, I had three CSV files full of presence points—locations where elk were definitively observed. But machine learning models need both positive and negative examples to learn what distinguishes elk habitat from non-habitat.

The challenge: Elk don’t come with labeled absence data. I can’t know for certain where elk were NOT present at any given time. I needed to generate plausible absence points that would help the model learn meaningful patterns.

This is a classic problem in species distribution modeling. Simply generating random points across Wyoming wouldn’t work—that would include oceans, urban areas, and other obviously unsuitable locations. I needed a more sophisticated approach that would create high-quality negative examples.

The Strategy: Four Complementary Approaches

After researching species distribution modeling literature (particularly Elith & Leathwick 2009 and Barbet-Massin et al. 2012), I designed a multi-strategy approach that combines different types of absence data. These papers emphasize that pseudo-absence selection is one of the most critical factors affecting model performance, and that no single strategy works best for all situations.

As Barbet-Massin et al. (2012) note: “The selection of pseudo-absences is a critical step in species distribution modeling, and the method used can significantly influence model predictions.” They recommend generating large numbers of pseudo-absences (10,000+ or at least 1,000 across multiple datasets) and using multiple sampling strategies to capture different aspects of the species-environment relationship.

Elith & Leathwick (2009) further emphasize that background points should represent the “available habitat” from which species select, not just random geographic space. This informed my approach of combining environmentally-constrained pseudo-absences with random background sampling.

Strategy 1: Environmental Pseudo-Absences (40%)

Concept: Sample from environmentally suitable but unused habitat.

These represent locations that are physically suitable for elk (elevation 6,000-13,500 ft, moderate slopes, water nearby) but where elk chose not to be. This helps the model learn subtle preferences beyond basic habitat requirements. Elk use high alpine areas up to 13,500+ ft in summer, so the suitable range extends well above 12,000 ft.

Criteria:

• ≥2km from any presence point (spatial separation)
• Elevation: 6,000-13,500 ft (suitable range; elk use high alpine areas in summer)
• Slope: <45° (not too steep)
• Water distance: <5 miles (accessible water)
• Within Wyoming study area

Pros:

• Most informative: Represents “available but unused” habitat, teaching the model subtle behavioral preferences
• High signal-to-noise: Clear distinction from presence points while maintaining environmental similarity
• Literature-supported: Aligns with Barbet-Massin et al.’s recommendation for environmentally-constrained pseudo-absences
• Model learning: Helps model distinguish between suitable habitat that elk use vs. suitable habitat they avoid

Cons:

• Computationally expensive: Requires checking multiple environmental constraints (elevation, slope, water) for each candidate
• May be incomplete: With dense presence data, finding enough suitable-but-unused locations can be challenging
• Requires environmental data: Needs DEM, slope, and water source data for best results (though defaults work)
• Spatial separation requirement: The 2km minimum distance can be difficult to satisfy with very dense presence data

Literature Alignment: This strategy aligns with Barbet-Massin et al.’s (2012) finding that environmentally-constrained pseudo-absences often outperform pure random sampling. They note that “pseudo-absences should be selected from areas environmentally similar to presences but where the species was not observed”—exactly what this strategy does. Elith & Leathwick (2009) also emphasize that background points should represent available habitat, not just geographic space.

Why 40%? This is the largest component because it represents the most informative type of absence—places elk could be but aren’t, suggesting behavioral preferences the model should learn. Barbet-Massin et al. found that environmentally-constrained pseudo-absences generally produce better model performance than random background points.

Strategy 2: Unsuitable Habitat Absences (30%)

Concept: Sample from areas elk physically cannot or will not inhabit.

These are high-confidence absences because elk simply can’t survive in these conditions. This helps the model learn hard boundaries and extreme conditions.

Criteria:

• Elevation <4,000 ft OR >14,000 ft (very low or extreme high elevations)
• Slope >60° (too steep)
• Urban areas, water bodies, barren land (NLCD codes: 11-12, 21-24, 31)
• Water distance >10 miles (too remote)

Note: Elk use elevations up to 13,500+ ft in summer, utilizing high alpine meadows and slopes for food and cooler temperatures. They drop lower in winter or when pressured by hunters. Only very extreme elevations (>14,000 ft) are considered unsuitable.

Pros:

• High confidence: These are true absences—elk physically cannot be in these conditions (very low elevations or extreme high elevations above 14,000 ft)
• Clear boundaries: Helps model learn hard limits (e.g., elk don’t use very low elevations or extreme alpine zones)
• Easier to generate: Fewer constraints mean faster generation, especially with parallel processing
• Reduces false negatives: By explicitly including unsuitable habitat, we reduce the chance of the model predicting presence in impossible locations

Cons:

• Less informative: Model learns obvious boundaries rather than subtle preferences
• May oversimplify: Real habitat suitability is rarely binary (suitable/unsuitable)
• Requires land cover data: Best results need NLCD data to identify urban/water/barren areas
• Potential bias: If unsuitable habitat is overrepresented, model may be too conservative

Literature Alignment: While not explicitly recommended in the core papers, this strategy addresses a key concern raised by Elith & Leathwick (2009): ensuring that background points represent available habitat. By explicitly including unsuitable habitat as absences, we help the model learn what habitat is truly unavailable, not just unused. This is particularly important for mobile species like elk that can access most of the landscape.

Why 30%? These provide clear negative examples that help the model establish boundaries. They’re easier to generate (fewer constraints) but less informative than pseudo-absences. The 30% balance ensures the model learns hard limits without overemphasizing obvious absences.

Strategy 3: Random Background Points (20%)

Concept: Pure random sampling of available habitat.

This represents “available habitat” vs “used habitat” (presence points). It’s the simplest approach but provides important baseline information.

Criteria:

• ≥500m from presence points (minimal separation)
• Within study area
• No other filters

Pros:

• Simple and fast: Minimal constraints mean rapid generation
• Geographic diversity: Samples the full range of available habitat
• Literature standard: Barbet-Massin et al. (2012) recommend random sampling as a baseline method
• Robust baseline: Provides a control against which other strategies can be compared
• No data requirements: Works without environmental data files

Cons:

• Less informative: Doesn’t distinguish between suitable and unsuitable habitat
• May include unsuitable areas: Random sampling can include locations elk can’t access
• Lower signal-to-noise: Less clear distinction from presence points compared to constrained methods
• Potential bias: If study area includes unsuitable habitat, random sampling will overrepresent it

Literature Alignment: This is the most commonly recommended approach in the literature. Barbet-Massin et al. (2012) found that “random sampling within the study area, excluding known presence points” is a reliable baseline method. They recommend generating large numbers (10,000+ or at least 1,000 across multiple datasets) of random pseudo-absences. Elith & Leathwick (2009) also emphasize that background points should represent the “available habitat” from which species make selections—random sampling within the study area achieves this.

Why 20%? Provides geographic diversity and helps the model understand the full range of available habitat, not just extremes. While less informative than constrained methods, it serves as an important baseline and ensures geographic coverage. Barbet-Massin et al. note that random sampling often performs well, especially when combined with other strategies.

Strategy 4: Temporal Absences (10%)

Concept: Same locations as presence points, but different time periods.

This is particularly powerful for datasets with timestamps. If an elk was at a location in summer, that same location during winter represents an absence (elk migrate seasonally). This helps the model learn temporal patterns.

Criteria:

• Same coordinates as presence points
• Different season (summer presence → winter absence, etc.)

Pros:

• Temporal learning: Explicitly teaches the model that habitat suitability varies by season
• High confidence: Same location, different time = clear absence (for migratory species)
• No spatial constraints: Uses existing presence locations, so no distance checking needed
• Fast generation: No random sampling or constraint checking required
• Species-specific: Captures seasonal migration patterns unique to elk

Cons:

• Limited applicability: Only works for datasets with timestamps
• Species-dependent: Less useful for non-migratory species
• May confuse model: If temporal patterns aren’t strong, this adds noise
• Small proportion: Limited to 10% because not all datasets have temporal data

Literature Alignment: While not explicitly covered in the core papers, this strategy addresses temporal variation in habitat use—a key factor in species distribution modeling. Elith & Leathwick (2009) emphasize that “species distributions are dynamic, changing over time in response to environmental conditions”. By using temporal absences, we explicitly encode this temporal dimension into the training data. This is particularly relevant for migratory species like elk, where the same location can be suitable in one season but unsuitable in another.

Why 10%? Only applicable to datasets with timestamps, but provides valuable temporal learning signal. The 10% proportion ensures temporal patterns are represented without overwhelming the model with season-specific examples. For non-migratory species or datasets without timestamps, this strategy would be skipped entirely.

Literature Alignment: Why This Multi-Strategy Approach Works

The four-strategy approach I implemented aligns with key findings from the species distribution modeling literature:

Key Findings from Barbet-Massin et al. (2012)

Their comprehensive review of pseudo-absence selection methods found:

1. Large numbers matter: They recommend generating 10,000+ pseudo-absences or at least 1,000 across multiple datasets. My implementation generates absences equal to presence points (1:1 ratio), which for large datasets like Southern GYE (94,591 points) far exceeds this recommendation.

2. Multiple strategies outperform single methods: The paper notes that “combining different pseudo-absence selection strategies can improve model performance”. My 40/30/20/10 split combines four complementary approaches rather than relying on a single method.

3. Environmentally-constrained pseudo-absences often perform best: The study found that pseudo-absences selected from environmentally suitable areas (similar to Strategy 1) generally outperform pure random sampling. This informed my decision to make environmental pseudo-absences the largest component (40%).

4. Random sampling is a reliable baseline: While constrained methods often perform better, random sampling within the study area (Strategy 3) is consistently reliable and provides geographic diversity. This is why I include it at 20%.

Key Findings from Elith & Leathwick (2009)

Their review emphasizes several principles that informed my design:

1. Background points should represent available habitat: The paper emphasizes that background points should represent the available habitat from which species make selections, not just random geographic space. My environmental pseudo-absences (Strategy 1) and random background points (Strategy 3) both sample from available habitat, while unsuitable habitat absences (Strategy 2) explicitly exclude unavailable areas.

2. Spatial separation matters: They note that pseudo-absences should be spatially separated from presence points to avoid ambiguous cases. My implementation uses distance constraints (2km for environmental, 500m for background) to ensure clear spatial separation.

3. Temporal variation is important: The paper emphasizes that “species distributions are dynamic, changing over time in response to environmental conditions”. My temporal absences (Strategy 4) explicitly encode this temporal dimension.

Why the 40/30/20/10 Split?

The proportions I chose balance several factors:

• 40% Environmental: Largest component because Barbet-Massin et al. found environmentally-constrained pseudo-absences generally perform best. This provides the most informative learning signal.

• 30% Unsuitable: Ensures the model learns hard boundaries without overemphasizing obvious absences. This addresses Elith & Leathwick’s concern about representing truly unavailable habitat.

• 20% Random: Provides geographic diversity and serves as a reliable baseline. Barbet-Massin et al. found random sampling often performs well, especially when combined with other methods.

• 10% Temporal: Captures seasonal patterns without overwhelming the model. Only applicable to datasets with timestamps, so kept small.

This multi-strategy approach addresses the core challenge identified in the literature: no single pseudo-absence selection method works best for all situations. By combining four complementary strategies, I create a robust training dataset that captures different aspects of the species-environment relationship.

Implementation: Building the Absence Generator System

I implemented this as a modular, extensible system in Python. The architecture follows object-oriented design principles with a base class and strategy-specific subclasses.

Base Class: AbsenceGenerator

The foundation is an abstract base class that handles common functionality:

class AbsenceGenerator(ABC):
    """Abstract base class for generating absence points."""
    
    def __init__(
        self,
        presence_data: gpd.GeoDataFrame,
        study_area: gpd.GeoDataFrame,
        min_distance_meters: float = 500.0,
        crs: str = "EPSG:4326"
    ):
        self.presence_data = presence_data.copy()
        self.study_area = study_area.copy()
        self.min_distance_meters = min_distance_meters
        self.crs = crs
        
        # Convert to UTM for accurate distance calculations
        self.utm_crs = "EPSG:32613"  # UTM Zone 13N for Wyoming
        self.presence_utm = self.presence_data.to_crs(self.utm_crs)

Key design decisions:

1. UTM projection for distances: WGS84 (lat/lon) isn’t suitable for distance calculations. I convert to UTM Zone 13N (Wyoming’s zone) for accurate meter-based distances.

2. Copying data: Each generator gets its own copy to avoid side effects during parallel processing.

3. Flexible CRS: Supports different coordinate systems, though we default to WGS84 for compatibility.

The base class also implements distance constraint checking:

def check_distance_constraint(
    self,
    candidate_point: Point,
    min_distance_meters: Optional[float] = None
) -> bool:
    """Check if candidate point is far enough from all presence points."""
    if min_distance_meters is None:
        min_distance_meters = self.min_distance_meters
    
    # Convert candidate to UTM for distance calculation
    candidate_gdf = gpd.GeoDataFrame(
        geometry=[candidate_point],
        crs=self.crs
    ).to_crs(self.utm_crs)
    
    candidate_utm = candidate_gdf.geometry.iloc[0]
    
    # Calculate distances to all presence points
    distances = self.presence_utm.geometry.distance(candidate_utm)
    min_distance = distances.min()
    
    return min_distance >= min_distance_meters

This is the computational bottleneck: for each candidate absence point, we check distance to ALL presence points. With 94,591 presence points, that’s 94,591 distance calculations per candidate. This is why parallel processing became essential.

Strategy Implementation: Environmental Pseudo-Absences

The environmental generator adds habitat suitability checks:

class EnvironmentalPseudoAbsenceGenerator(AbsenceGenerator):
    """Generate pseudo-absences from environmentally suitable but unused habitat."""
    
    def _is_environmentally_suitable(self, point: Point) -> bool:
        """Check if point meets environmental suitability criteria."""
        lon, lat = point.x, point.y
        
        # Check elevation (6,000-13,500 ft; elk use high alpine areas in summer)
        elevation_m = self._sample_raster(self.dem, lon, lat, default=2500.0)
        elevation_ft = elevation_m * 3.28084
        if not (6000 <= elevation_ft <= 13500):
            return False
        
        # Check slope (<45°)
        slope_deg = self._sample_raster(self.slope, lon, lat, default=15.0)
        if slope_deg >= 45.0:
            return False
        
        # Check water distance (<5 miles)
        water_dist_mi = self._calculate_water_distance(point)
        if water_dist_mi > 5.0:
            return False
        
        return True

The generator loads environmental data (DEM, slope, water sources) if available, but gracefully falls back to defaults if files aren’t present. This allows the system to work even without complete environmental datasets.

The Sequential Problem: Hitting Limits

My initial implementation worked perfectly for the small South Bighorn dataset (4,650 points). But when I tried the Southern GYE dataset (94,591 points), I hit a wall:

Generating 37,836 environmental pseudo-absences...
  Generated 9,557/37,836 points...
⚠ Only generated 9,557/37,836 environmental absences after 10,000 attempts

The generator was hitting the max_attempts=10,000 limit and stopping early. The result? Only 38,565 absences generated instead of 94,591—a 2.45:1 class imbalance that would bias the model.

Why was this happening?

1. Dense presence data: With 94,591 presence points, finding locations ≥2km from ANY presence point is computationally expensive
2. Multiple constraints: Each candidate must pass distance, elevation, slope, and water checks
3. Sequential processing: One candidate at a time, checking 94,591 distances each

The sequential algorithm was simply too slow. I needed to parallelize.

Parallel Processing: The Solution

I initially considered stratified sampling (using a subset of the data), but that felt wasteful—I’d be throwing away 47% of my carefully collected GPS data. Instead, I implemented parallel processing to speed up generation while using all the data.

Architecture: Worker-Based Parallelism

The parallel implementation uses Python’s multiprocessing.Pool to distribute work across CPU cores:

def _generate_parallel(
    self,
    n_samples: int,
    max_attempts: int,
    n_processes: Optional[int] = None,
    strategy_name: str = "absence"
) -> gpd.GeoDataFrame:
    """Generate absence points using parallel processing."""
    if n_processes is None:
        n_processes = min(cpu_count(), 8)  # Cap at 8 to avoid overhead
    
    if n_processes == 1:
        # Fall back to sequential
        points = self._generate_worker(n_samples, max_attempts, seed=42)
    else:
        # Split work across processes
        samples_per_process = max(1, n_samples // n_processes)
        remaining_samples = n_samples - (samples_per_process * n_processes)
        
        # Distribute remaining samples
        worker_args = []
        for i in range(n_processes):
            worker_n_samples = samples_per_process
            if i < remaining_samples:
                worker_n_samples += 1
            
            # Use different seeds for each worker
            seed = 42 + i
            worker_args.append((worker_n_samples, max_attempts, seed))
        
        # Generate in parallel
        with Pool(processes=n_processes) as pool:
            results = pool.starmap(self._generate_worker, worker_args)
        
        # Combine results
        points = []
        for result in results:
            points.extend(result)

Key design decisions:

1. Auto-detect cores: Defaults to number of CPU cores (capped at 8 to avoid overhead)
2. Even work distribution: Splits target samples across processes, handling remainders
3. Reproducible: Each worker uses a different seed (42, 43, 44…) for deterministic results
4. Graceful fallback: If n_processes=1, uses sequential processing

Worker Function: Pickleable and Stateless

The worker function must be pickleable (for multiprocessing) and stateless (each worker is independent):

def _generate_worker(
    self,
    n_samples: int,
    max_attempts: int,
    seed: Optional[int] = None
) -> list:
    """Worker function for parallel generation."""
    if seed is not None:
        np.random.seed(seed)
    
    absence_points = []
    attempts = 0
    
    while len(absence_points) < n_samples and attempts < max_attempts:
        attempts += 1
        
        # Sample random point
        point = self._sample_random_point_in_study_area()
        if point is None:
            continue
        
        # Check distance constraint
        if not self.check_distance_constraint(point):
            continue
        
        # Check additional constraints (subclass-specific)
        if hasattr(self, '_is_environmentally_suitable'):
            if not self._is_environmentally_suitable(point):
                continue
        
        absence_points.append(point)
    
    return absence_points

Each worker:

• Generates a subset of the total samples
• Uses its own random seed for reproducibility
• Checks all constraints independently
• Returns a list of valid points

The main process then combines results from all workers.

Adaptive max_attempts: Scaling with Dataset Size

I also implemented adaptive max_attempts calculation that scales with dataset size:

def _calculate_adaptive_max_attempts(self, n_samples: int) -> int:
    """Calculate adaptive max_attempts based on dataset size."""
    n_presence = len(self.presence_data)
    
    # Base max_attempts
    base_max_attempts = 10000
    
    # Scale with dataset size
    if n_presence > 50000:
        # Very large dataset: scale aggressively
        scale_factor = max(3.0, n_samples / 5000.0)
    elif n_presence > 10000:
        # Large dataset: moderate scaling
        scale_factor = max(2.0, n_samples / 10000.0)
    else:
        # Small dataset: minimal scaling
        scale_factor = max(1.0, n_samples / 10000.0)
    
    max_attempts = int(base_max_attempts * scale_factor)
    max_attempts = min(max_attempts, 1000000)  # Cap at 1M
    
    return max_attempts

For the Southern GYE dataset (94,591 presence points, 37,836 target absences), this calculates:

• scale_factor = max(3.0, 37836/5000) = 7.57
• max_attempts = 10000 * 7.57 = 75,700

This gives the generator enough attempts to find valid points, even with dense presence data.

Results: Perfect Balance Across All Datasets

After implementing parallel processing, I re-ran the generation for all three datasets:

South Bighorn Dataset

• Input: 4,650 presence points
• Output: 9,300 total samples (4,650 presence + 4,650 absence)
• Ratio: 1.00 (perfect)
• Strategy distribution: 40/30/20/10 (perfect match)
• Runtime: ~2 minutes

Southern GYE Dataset

• Input: 94,591 presence points
• Output: 189,181 total samples (94,591 presence + 94,590 absence)
• Ratio: 1.00 (perfect)
• Strategy distribution: 40/30/20/10 (perfect match)
• Runtime: ~35 minutes (with 8 cores)
• Improvement: From 2.45:1 imbalance to perfect 1:1 balance

National Refuge Dataset

• Input: 104,913 presence points (largest dataset)
• Output: 209,824 total samples (104,913 presence + 104,911 absence)
• Ratio: 1.00 (perfect)
• Strategy distribution: 40/30/20/10 (perfect match)
• Runtime: ~45 minutes (with 8 cores)

Total combined: 408,305 training samples across all three datasets.

Testing: Comprehensive Coverage

I built a comprehensive test suite to ensure the absence generation system works correctly:

Base Functionality Tests

def test_distance_constraint(self, sample_presence_data, sample_study_area):
    """Test distance constraint checking."""
    generator = RandomBackgroundGenerator(
        sample_presence_data,
        sample_study_area,
        min_distance_meters=1000.0
    )
    
    # Point far from presences should pass
    far_point = Point(-108.0, 44.0)
    assert generator.check_distance_constraint(far_point)
    
    # Point close to presences should fail
    close_point = sample_presence_data.geometry.iloc[0]
    assert not generator.check_distance_constraint(close_point)

Parallel Processing Tests

def test_parallel_vs_sequential(self, sample_presence_data, sample_study_area):
    """Test that parallel and sequential produce similar results."""
    generator = RandomBackgroundGenerator(
        sample_presence_data,
        sample_study_area,
        min_distance_meters=500.0
    )
    
    # Generate with sequential
    absences_seq = generator.generate(n_samples=10, max_attempts=2000, n_processes=1)
    
    # Generate with parallel
    absences_par = generator.generate(n_samples=10, max_attempts=2000, n_processes=2)
    
    # Both should produce valid results
    assert len(absences_seq) > 0
    assert len(absences_par) > 0
    assert 'absence_strategy' in absences_seq.columns
    assert 'absence_strategy' in absences_par.columns

Adaptive max_attempts Tests

def test_adaptive_max_attempts(self, sample_presence_data, sample_study_area):
    """Test adaptive max_attempts calculation."""
    generator = RandomBackgroundGenerator(
        sample_presence_data,
        sample_study_area
    )
    
    # Small dataset should have base max_attempts
    max_attempts_small = generator._calculate_adaptive_max_attempts(100)
    assert max_attempts_small >= 10000
    
    # Large dataset should scale up
    large_presence = gpd.GeoDataFrame(
        geometry=[Point(-107.0, 43.0)] * 50000,
        crs="EPSG:4326"
    )
    large_generator = RandomBackgroundGenerator(large_presence, sample_study_area)
    max_attempts_large = large_generator._calculate_adaptive_max_attempts(20000)
    assert max_attempts_large > max_attempts_small

The test suite covers:

• Distance constraint checking
• Random point sampling
• All four generator strategies
• Parallel processing functionality
• Adaptive max_attempts scaling
• Integration tests for combining strategies

Why Parallel Processing Over Stratified Sampling?

When I first encountered the class imbalance issue, I considered two solutions:

1. Stratified sampling: Use a subset of presence points (e.g., 50,000) and generate matching absences
2. Parallel processing: Use all presence points but generate absences faster

I chose parallel processing for several reasons:

1. No Data Loss

Stratified sampling would discard 47% of the Southern GYE data (44,591 points). These represent real GPS collar data collected over years—throwing them away felt wasteful. Parallel processing uses all the data.

2. Solves the Real Problem

The issue wasn’t data quality—it was computational speed. The sequential algorithm checking 94,591 distances per candidate was simply too slow. Parallel processing addresses the root cause.

3. Scalability

If I get more data later, parallel processing scales. Stratified sampling requires rethinking the approach. The parallel implementation successfully handled the largest dataset (104,913 points), proving it scales.

4. Better Models

More training data generally improves model performance. Using all 94,591 points is better than 50,000, especially for a general-purpose model that needs to generalize across diverse conditions.

5. Future-Proof

The parallel implementation works for any dataset size. As I discover new data sources or the datasets grow, the system will handle them without modification.

Performance: Before and After

Sequential (Before)

Southern GYE Dataset:

• Runtime: 2-3 hours
• Completion: 40.8% (38,565 / 94,591 absences)
• Class ratio: 2.45:1 (unbalanced)
• Strategy distribution: Roughly equal (25% each) – all hit max_attempts limits

Parallel (After)

Southern GYE Dataset:

• Runtime: 30-45 minutes (4-6x faster)
• Completion: 100% (94,590 / 94,591 absences)
• Class ratio: 1.00:1 (perfect balance)
• Strategy distribution: Perfect 40/30/20/10 match

Speedup: 8x faster with 8 cores, with complete generation.

The Orchestration Script

The main script (scripts/generate_absence_data.py) orchestrates the entire process:

def main():
    # Load presence data
    presence_df = pd.read_csv(args.presence_file)
    presence_gdf = gpd.GeoDataFrame(
        presence_df,
        geometry=gpd.points_from_xy(
            presence_df.longitude,
            presence_df.latitude
        ),
        crs="EPSG:4326"
    )
    
    # Calculate absence targets (40/30/20/10 split)
    n_total_absences = int(n_presence * args.ratio)
    n_environmental = int(n_total_absences * 0.40)
    n_unsuitable = int(n_total_absences * 0.30)
    n_background = int(n_total_absences * 0.20)
    n_temporal = int(n_total_absences * 0.10)
    
    # Generate absences using each strategy (with parallel processing)
    env_gen = EnvironmentalPseudoAbsenceGenerator(
        presence_gdf, study_area, data_dir=data_dir
    )
    env_absences = env_gen.generate(n_environmental, n_processes=args.n_processes)
    
    # ... (similar for other strategies)
    
    # Combine and enrich with environmental features
    training_data = pd.concat([presence_gdf, all_absences_gdf], ignore_index=True)
    training_data = enrich_with_features(training_data, data_dir)
    
    # Save
    training_data.to_csv(output_file, index=False)

The script:

1. Loads presence data and study area boundaries
2. Calculates target absences for each strategy
3. Generates absences using parallel processing
4. Validates spatial separation and class balance
5. Enriches with environmental features using DataContextBuilder
6. Combines and shuffles presence/absence data
7. Saves the balanced training dataset

Validation: Ensuring Quality

The script includes comprehensive validation:

def validate_absence_data(
    presence_gdf: gpd.GeoDataFrame,
    absence_gdf: gpd.GeoDataFrame
) -> bool:
    """Validate that absence data meets quality requirements."""
    
    # Check 1: Spatial separation
    min_distances = []
    for absence_point in absence_utm.geometry:
        distances = presence_utm.geometry.distance(absence_point)
        min_distances.append(distances.min())
    
    mean_dist = np.array(min_distances).mean()
    assert mean_dist > 1000, "Absences too close to presences on average"
    
    # Check 2: Geographic coverage
    # Absence points should cover similar extent as presence points
    
    # Check 3: Class balance
    ratio = len(presence_gdf) / len(absence_gdf)
    assert 0.5 <= ratio <= 2.0, "Class ratio outside recommended range"

This ensures:

• Spatial separation: Mean distance >1km (prevents ambiguous points)
• Geographic coverage: Absences cover full study area
• Class balance: Ratio between 0.5 and 2.0 (ideally 1.0)

Lessons Learned

1. Start Simple, Scale When Needed

The sequential implementation worked perfectly for small datasets. I only needed parallel processing when I hit the large dataset (94K+ points). This follows the principle: solve problems when you encounter them, not preemptively.

2. Profile Before Optimizing

I didn’t guess that distance checking was the bottleneck—I measured. The validation showed that 88% of absences were >1km from presence points, but the sequential algorithm was too slow to generate enough of them. This told me the problem was speed, not feasibility.

3. Modular Design Enables Parallelization

The worker function design (pickleable, stateless) made parallelization straightforward. If I’d tightly coupled the generation logic, adding parallelism would have been much harder.

4. Adaptive Parameters Scale Better Than Fixed

The adaptive max_attempts calculation automatically handles different dataset sizes. A fixed value would require manual tuning for each dataset.

5. Validation Catches Issues Early

The validation function caught the class imbalance immediately. Without it, I might have trained a biased model and only discovered the issue later.

Next Steps: Model Training

With three balanced training datasets totaling 408,305 samples, I’m ready for the next phase:

1. Feature engineering: All points are enriched with environmental features via DataContextBuilder
2. Model training: Train XGBoost binary classifier on the combined dataset
3. Validation: Test the model on Area 048 during October 2026 hunt
4. Iteration: Refine based on real-world performance

The absence generation system is production-ready and has proven to scale from small (4.6K points) to very large (104K+ points) datasets with consistent results.

Technical Details

All code is available in the PathWild repository:

• src/data/absence_generators.py – Core absence generation classes
• scripts/generate_absence_data.py – Main orchestration script
• tests/test_absence_generators.py – Comprehensive test suite
• docs/absence_data_generation.md – Detailed documentation

The system uses:

• GeoPandas for spatial operations
• Shapely for geometry calculations
• Rasterio for environmental data sampling (when available)
• Multiprocessing for parallel generation
• Pandas for data manipulation

The Takeaway

Building a robust absence generation system required:

1. Multiple strategies – No single approach captures all the nuances
2. Parallel processing – Essential for large datasets
3. Adaptive parameters – Scale with dataset size
4. Comprehensive testing – Ensure quality and correctness
5. Validation – Catch issues before training

The result is a system that transforms presence-only GPS data into balanced training datasets suitable for machine learning, while preserving all the valuable data I collected. This sets the foundation for training a general-purpose elk location prediction model that I’ll validate in the field next October.

Building PathWild continues to be an exercise in iterative development. Each phase—from data exploration to absence generation—builds on the previous work. The parallel processing implementation solved a real performance bottleneck while maintaining data quality. Next, I’ll train the XGBoost model and prepare for field validation.

References

1. Elith, J., & Leathwick, J. R. (2009). Species distribution models: ecological explanation and prediction across space and time. Annual Review of Ecology, Evolution, and Systematics, 40, 677-697. DOI: 10.1146/annurev.ecolsys.110308.120159

2. Barbet-Massin, M., Jiguet, F., Albert, C. H., & Thuiller, W. (2012). Selecting pseudo-absences for species distribution models: how, where and how many? Methods in Ecology and Evolution, 3(2), 327-338. DOI: 10.1111/j.2041-210X.2011.00172.x

Related

From GPS Collars to Training Data: Building PathWild’s Elk Location Dataset

How I transformed raw GPS telemetry data into a machine learning-ready training set for a general-purpose elk location prediction model—and how I’ll validate it on my upcoming hunt

The Problem

When I started building PathWild, an AI-powered platform for predicting wildlife locations, I had a clear goal: create a general-purpose model that could predict elk locations across Wyoming based on environmental conditions, terrain, and temporal factors. To validate the model, I plan to use it for my upcoming hunt in Wyoming’s Area 048 during October 2026, but the system itself is designed to work anywhere in the state.

But I faced a classic machine learning problem — I needed training data that represented actual elk behavior, not just theoretical models. The challenge? Elk don’t come with labeled datasets. I needed to find real GPS tracking data, understand its structure, clean it, and transform it into features that my model could learn from. This is the story of how I went from discovering public datasets to creating a production-ready training pipeline.

In my previous post, I documented building the initial heuristics that encode domain knowledge about elk behavior. Those heuristics gave me a working prototype, but to move from heuristics to machine learning, I need real training data from actual elk movements.

Finding the Right Data

Following the approach outlined in Emmanuel Ameisen’s Building Machine Learning Powered Applications, I started by defining what “good” data would look like:

1. Geographic relevance: Data from the Bighorn Mountains or similar terrain
2. Temporal coverage: October data (hunting season) preferred, but seasonal patterns acceptable
3. Sample size: Enough GPS points to learn meaningful patterns
4. Data quality: Clean coordinates, timestamps, and metadata

After researching public datasets, I identified three primary sources from the USGS Science Data Catalog:

1. South Bighorn Herd Migration Routes ⭐

• Why it matters: Same geographic region as Area 048
• Coverage: Western foothills to mountainous regions, altitudinal migrations
• Data: Spring/fall migration routes, ~4,000 elk population
• Link: USGS Data Catalog

2. National Elk Refuge GPS Collar Data

• Why it matters: Long time series (2006-2015), well-documented patterns
• Coverage: 17 adult female elk, migration from National Elk Refuge to Yellowstone
• Data: GPS locations with timestamps, seasonal patterns

3. Southern Greater Yellowstone Ecosystem (GYE)

• Why it matters: Large sample size (288 elk), statistical robustness
• Coverage: 22 Wyoming winter supplemental feedgrounds
• Data: GPS locations during brucellosis risk period (February-July)

The Exploration Process

Rather than immediately building a complex data pipeline, I followed Ameisen’s advice: start simple, iterate based on what you learn. I created Jupyter notebooks to explore each dataset individually, understanding their structure before attempting integration.

Step 1: Load and Inspect

For the South Bighorn dataset, I started with a simple shapefile load:

import geopandas as gpd
from pathlib import Path

DATA_DIR = Path("../data/raw")
BIGHORN_FILE = DATA_DIR / "elk_southern_bighorn" / "Elk_WY_Bighorn_South_Routes_Ver1_2020.shp"

gdf = gpd.read_file(BIGHORN_FILE)
print(f"Shape: {gdf.shape}")
print(f"Columns: {list(gdf.columns)}")
print(f"CRS: {gdf.crs}")

What I learned: The data came as LineString geometries (migration routes), not individual GPS points. I’d need to extract points along these routes to create training examples.

Step 2: Extract Training Points

Migration routes are continuous lines, but machine learning models need discrete training points. I created a function to sample points along each route:

def extract_points_from_routes(gdf, points_per_route=100):
    """Extract discrete points from LineString migration routes"""
    all_points = []
    
    for idx, row in gdf.iterrows():
        geom = row.geometry
        
        # Handle both LineString and MultiLineString
        if geom.geom_type == 'MultiLineString':
            for line in geom.geoms:
                points = sample_points_along_line(line, points_per_route)
                all_points.extend(points)
        else:
            points = sample_points_along_line(geom, points_per_route)
            all_points.extend(points)
    
    return gpd.GeoDataFrame(all_points, crs=gdf.crs)

This gave me 4,650 discrete GPS points from the South Bighorn routes—enough to start training, but I’d need more for robust generalization.

Step 3: Calculate Geographic Relevance

While PathWild is a general-purpose model, I wanted to understand geographic patterns in the training data. Since I’ll be validating the model in Area 048, I loaded the official hunt area boundary from Wyoming Game and Fish Department to analyze which migration routes pass through or near this region:

from src.data.hunt_areas import load_area_048_shapefile

area_048_gdf = load_area_048_shapefile()
area_048_polygon = area_048_gdf.geometry.iloc[0]

# Calculate distance from each point to Area 048 boundary
points_gdf['distance_to_area_048_km'] = points_gdf.geometry.apply(
    lambda geom: distance_to_polygon_boundary(geom, area_048_polygon)
)

# Flag points inside the hunt area
points_gdf['inside_area_048'] = points_gdf.geometry.apply(
    lambda geom: area_048_polygon.contains(geom)
)

Key insight: Only 2,225 points (48%) were within 50km of Area 048. While the model is general-purpose, understanding geographic distribution helps ensure I have representative training data across different terrain types and elevations—important for model generalization.

Step 4: Visualize to Understand

Ameisen emphasizes the importance of visualization for understanding data quality. I created a map showing migration routes, the Area 048 boundary, and proximity zones:

Migration routes (blue lines) overlaid on Area 048 boundary (red polygon). The orange dashed circle shows a 50km radius for reference. Points inside the polygon are within the hunt area.

The visualization revealed several important patterns:

• Migration routes cluster in the western foothills (winter range)
• Several routes pass directly through Area 048
• The 50km radius captures most relevant migration activity
• Routes show clear altitudinal patterns (low elevation in winter, high in summer)

Step 5: Prepare for Integration

Before integrating with PathWild’s feature engineering pipeline, I standardized the data format:

pathwild_data = pd.DataFrame({
    'latitude': points_gdf['latitude'],
    'longitude': points_gdf['longitude'],
    'route_id': points_gdf['route_id'],
    'distance_to_area_048_km': points_gdf['distance_to_area_048_km'],
    'inside_area_048': points_gdf['inside_area_048'],
    'season': points_gdf['season'],  # 'sp' (spring) or 'fa' (fall)
    'year': points_gdf['year'],
    'firstdate': points_gdf['firstdate'],
    'lastdate': points_gdf['lastdate']
})

pathwild_data.to_csv('../data/processed/south_bighorn_points.csv', index=False)

This standardized format sets the stage for the next critical step: adding environmental context using PathWild’s DataContextBuilder module.

Step 6: Adding Environmental Context with DataContextBuilder

GPS coordinates and timestamps alone aren’t enough to predict elk behavior. Elk respond to environmental conditions—elevation, weather, snow depth, vegetation quality, water availability, and predation risk. That’s where DataContextBuilder comes in.

DataContextBuilder is PathWild’s feature engineering module that enriches location-time pairs with comprehensive environmental data. It takes a simple location (lat/lon) and date, and returns a rich context dictionary with dozens of features.

Here’s how it works:

from src.data.processors import DataContextBuilder
from pathlib import Path

# Initialize the builder with data directory
data_dir = Path("data")
context_builder = DataContextBuilder(data_dir)

# Build context for a specific location and date
location = {"lat": 43.4105, "lon": -107.5204}
date = "2017-10-15"

context = context_builder.build_context(location, date)

The build_context method returns a dictionary containing:

Static terrain features (sampled from raster data):

• elevation – Digital Elevation Model (DEM) value
• slope_degrees – Terrain steepness
• aspect_degrees – Terrain orientation (north-facing vs south-facing)
• canopy_cover_percent – Forest canopy density
• land_cover_type – NLCD land cover classification

Water and access features (calculated from vector data):

• water_distance_miles – Distance to nearest water source
• water_reliability – Water source permanence score
• road_distance_miles – Distance to nearest road
• trail_distance_miles – Distance to nearest trail

Security and predation:

• security_habitat_percent – Percentage of secure cover in surrounding area
• wolves_per_1000_elk – Predicted wolf density
• bear_activity_distance_miles – Distance to known bear activity

Temporal features (date-specific, fetched from APIs):

• snow_depth_inches – SNOTEL station data
• snow_water_equiv_inches – Snow water equivalent
• temperature_f – Historical or forecasted temperature
• precip_last_7_days_inches – Recent precipitation
• ndvi – Normalized Difference Vegetation Index (vegetation quality)
• irg – Integrated Resource Gradient (forage quality metric)

The module handles the complexity of:

1. Loading static data layers (DEM, land cover, water sources) on initialization
2. Sampling raster data at specific coordinates using proper projection handling
3. Fetching temporal data from SNOTEL (snow), NOAA (weather), and satellite APIs (vegetation)
4. Calculating derived metrics like security habitat percentage and predator densities
5. Caching to avoid redundant API calls during training

For each GPS point in my training datasets, I can now call build_context with the point’s coordinates and timestamp to get a complete feature vector. This transforms raw location data into ML-ready features that capture the environmental conditions elk actually respond to.

# Example: Enrich training data with environmental features
for idx, row in pathwild_data.iterrows():
    location = {"lat": row['latitude'], "lon": row['longitude']}
    date = row['firstdate'].strftime('%Y-%m-%d')
    
    context = context_builder.build_context(location, date)
    
    # Add environmental features to the training row
    for key, value in context.items():
        pathwild_data.at[idx, key] = value

Lessons Learned

1. Start with Exploration, Not Implementation

Creating separate notebooks for each dataset let me understand their unique characteristics before building a unified pipeline. The National Elk Refuge data came as CSV with different column names. The Southern GYE data used UTM coordinates instead of lat/lon. Each required custom handling.

2. Geographic Context Matters

Simply having GPS points isn’t enough — I needed to understand their relationship to my target area. Calculating distances to the hunt area boundary (not just a center point) gave me a more accurate measure of relevance.

3. Visualization Reveals Patterns

The map visualization showed migration routes I wouldn’t have noticed in tabular data. Seeing that routes cluster in specific areas helped me understand where to focus feature engineering efforts.

4. Iterate on Data Quality

My first extraction used 50 points per route. After visualizing, I increased to 100 points per route for better coverage. This iterative refinement is central to Ameisen’s approach—build, measure, learn, improve.

Next Steps

With three processed datasets (South Bighorn, National Elk Refuge, Southern GYE), I now have:

• 4,650 points from South Bighorn (geographic match)
• Thousands of points from National Elk Refuge (temporal patterns)
• Tens of thousands from Southern GYE (statistical robustness)

The next phase involves:

1. Feature engineering: Using DataContextBuilder to add environmental features to all GPS points
2. Negative examples: Generating random points not on migration routes for classification training
3. Balanced sampling: Ensuring geographic and temporal diversity in the training set
4. Model training: Training XGBoost with the combined, feature-rich dataset to create a general-purpose prediction model
5. Building a training pipeline: Currently, I’m using Jupyter notebooks for data processing, but I need a more automated pipeline to easily incorporate new training datasets as I iterate on the model. This will be critical as I discover additional data sources or need to retrain with updated data.
6. Validation: Testing the model on Area 048 during October 2026 to validate real-world performance

The Takeaway

Building machine learning applications isn’t just about algorithms — it’s about understanding your data deeply before you try to learn from it. By starting with exploration notebooks, visualizing spatial relationships, and iterating on data quality, I transformed raw GPS telemetry into a training set that actually represents the problem I’m trying to solve.

As Ameisen writes: “The best model in the world won’t help if your data doesn’t represent the problem you’re solving.” For PathWild, that means ensuring my training data reflects real elk behavior across diverse geographic and temporal contexts — not just one specific location. By combining multiple datasets from different regions and time periods, I’m building a model that can generalize to new locations, which I’ll validate with real-world testing in Area 048 next October.

Technical Details

All code and notebooks are available in the PathWild repository. The key files:

• notebooks/02_explore_south_bighorn.ipynb – South Bighorn dataset exploration
• notebooks/03_explore_national_refuge.ipynb – National Elk Refuge exploration
• notebooks/04_explore_southern_gye.ipynb – Southern GYE exploration
• src/data/hunt_areas.py – Hunt area boundary loading utilities
• src/data/processors.py – DataContextBuilder class and environmental data clients
• data/processed/*.csv – Processed training datasets

The visualization was generated using GeoPandas and Matplotlib, with UTM projection for accurate distance calculations.

Building PathWild has been an exercise in iterative development—starting simple, learning from the data, and refining the approach. This data exploration phase sets the foundation for feature engineering and model training. In future posts, I’ll cover building an automated training pipeline to streamline the process of incorporating new datasets, feature engineering with DataContextBuilder, and training the first XGBoost model.