analytics.py

"""
analytics.py
============
Statistical analysis, nucleotide composition, GC-content profiling, and comparative genomics.
"""

from collections import Counter
from typing import Dict, List, Tuple, Union
import pandas as pd


def calculate_gc_content(seq: str) -> float:
    """
    Calculate the GC-content percentage of a nucleotide sequence.
    
    Formula: (Count(G) + Count(C)) / Total_Length * 100
    
    Parameters
    ----------
    seq : str
        DNA or RNA nucleotide sequence.
        
    Returns
    -------
    float
        GC percentage (0.0 to 100.0). Returns 0.0 for empty sequences.
    """
    if not seq:
        return 0.0
    seq_upper = seq.upper()
    g_count = seq_upper.count("G")
    c_count = seq_upper.count("C")
    return ((g_count + c_count) / len(seq_upper)) * 100.0


def calculate_base_composition(seq: str) -> Dict[str, Union[int, float]]:
    """
    Calculate absolute counts and relative frequencies of all bases.
    
    Parameters
    ----------
    seq : str
        DNA or RNA nucleotide sequence.
        
    Returns
    -------
    Dict[str, Union[int, float]]
        Dictionary containing counts, percentages, and total sequence length.
    """
    total = len(seq)
    if total == 0:
        return {"total_length": 0, "gc_content": 0.0}

    seq_upper = seq.upper()
    counts = Counter(seq_upper)
    
    results: Dict[str, Union[int, float]] = {
        "total_length": total,
        "count_A": counts.get("A", 0),
        "count_T": counts.get("T", 0),
        "count_C": counts.get("C", 0),
        "count_G": counts.get("G", 0),
        "count_U": counts.get("U", 0),
        "count_other": total - (
            counts.get("A", 0) + counts.get("T", 0) + 
            counts.get("C", 0) + counts.get("G", 0) + 
            counts.get("U", 0)
        ),
        "freq_A_%": round((counts.get("A", 0) / total) * 100, 2),
        "freq_T_%": round((counts.get("T", 0) / total) * 100, 2),
        "freq_C_%": round((counts.get("C", 0) / total) * 100, 2),
        "freq_G_%": round((counts.get("G", 0) / total) * 100, 2),
        "freq_U_%": round((counts.get("U", 0) / total) * 100, 2),
        "gc_content_%": round(calculate_gc_content(seq_upper), 2),
    }
    return results


def sliding_window_gc(
    seq: str,
    window_size: int = 500,
    step: int = 100
) -> Tuple[List[int], List[float]]:
    """
    Compute GC-content along a sequence using a sliding window.
    
    Parameters
    ----------
    seq : str
        Genomic sequence.
    window_size : int, default 500
        Length of the sliding window in base pairs.
    step : int, default 100
        Step size (stride) to advance the window.
        
    Returns
    -------
    Tuple[List[int], List[float]]
        (positions, gc_values) where positions are window center coordinates.
    """
    seq_upper = seq.upper()
    seq_len = len(seq_upper)
    if seq_len < window_size:
        return [seq_len // 2], [calculate_gc_content(seq_upper)]

    positions: List[int] = []
    gc_values: List[float] = []

    for start in range(0, seq_len - window_size + 1, step):
        subseq = seq_upper[start:start + window_size]
        center = start + (window_size // 2)
        positions.append(center)
        gc_values.append(calculate_gc_content(subseq))

    return positions, gc_values


def compare_sequences(records: Dict[str, str]) -> pd.DataFrame:
    """
    Produce a summary DataFrame comparing multiple sequences.
    
    Parameters
    ----------
    records : Dict[str, str]
        Mapping of sequence identifier -> sequence string.
        
    Returns
    -------
    pd.DataFrame
        Pandas DataFrame of genomic metrics for each sequence.
    """
    rows = []
    for name, seq in records.items():
        comp = calculate_base_composition(seq)
        comp["identifier"] = name
        rows.append(comp)
        
    df = pd.DataFrame(rows)
    cols = ["identifier", "total_length", "gc_content_%", "freq_A_%", "freq_T_%", "freq_C_%", "freq_G_%", "count_other"]
    present_cols = [c for c in cols if c in df.columns]
    return df[present_cols]