The Full Moon Myth Takes the Field – Marcos Huerta

import itertools


def make_comeback_expr(inning=6, lead=2):
    return (
        pl.col(f'comeback_after_{inning}').eq(-1)
        & (pl.col(f'score_after_{inning}_h') - pl.col(f'score_after_{inning}_v'))
        .abs()
        .ge(lead)
    ).alias(f'comeback_of_{lead}_after_{inning}')


teamstats = pl.read_csv('teamstats.csv', infer_schema_length=100000, ignore_errors=True)
gid_comeback_lookup = (
    teamstats
    .filter(pl.col('stattype').eq('value'))
    .with_columns(score_after_5=pl.sum_horizontal([f'inn{n}' for n in range(1, 6)]))
    .with_columns(score_after_6=pl.sum_horizontal([f'inn{n}' for n in range(1, 7)]))
    .with_columns(score_after_7=pl.sum_horizontal([f'inn{n}' for n in range(1, 8)]))
    .with_columns(score_after_8=pl.sum_horizontal([f'inn{n}' for n in range(1, 9)]))
    .select([
        'gid',
        'vishome',
        'win',
        'score_after_5',
        'score_after_6',
        'score_after_7',
        'score_after_8',
    ])
    .pivot(
        on='vishome',
        values=['score_after_5', 'score_after_6', 'score_after_7', 'score_after_8', 'win'],
        index='gid',
    )
    .with_columns(*[
        (
            (pl.col(f'score_after_{n}_h') - pl.col(f'score_after_{n}_v')).sign()
            * (pl.col('win_h') - pl.col('win_v')).sign()
        ).alias(f'comeback_after_{n}')
        for n in range(5, 9)
    ])
    .with_columns(*[
        make_comeback_expr(i, l) for i, l in itertools.product([5, 6, 7, 8], [2, 3])
    ])
)

homer_data = pl.read_csv('savant_data.csv').with_columns(
    date=pl.col('game_date').str.to_date('%Y-%m-%d')
)
moon_data = make_moon()
merged = (
    homer_data
    .join(
        moon_data,
        on='date',
    )
    .with_columns(*[
        pl.col('hit_distance_sc').ge(d).alias(f'homer_over_{d}') for d in range(380, 470, 10)
    ])
    .with_row_index('gid')
)

homers_by_distance = pd.concat([
    make_phase_bin_data(merged, phase_range=[99], stat=stat)[['ratio', 'ratio_error']].assign(
        stat_name=stat, pct_increase=pd.col('ratio') - 1
    )
    for stat in [
        'homer_over_380',
        'homer_over_390',
        'homer_over_400',
        'homer_over_410',
        'homer_over_420',
        'homer_over_430',
        'homer_over_440',
        'homer_over_450',
        'homer_over_460',
        #    'homer_over_470',
    ]
])

fig = px.bar(
    pl.from_pandas(homers_by_distance).with_columns(
        pl
        .col('stat_name')
        .str.replace_all('_', ' ', literal=True)
        .str.replace_all('homer', 'Homers')
        .str.to_titlecase()
        .str.replace_all(r'Of (\d)', r'Of $1+')
    ),
    y='stat_name',
    x='pct_increase',
    error_x='ratio_error',
    labels={
        'pct_increase': 'Percent Increase during Full Moon',
        'stat_name': 'Home run distance',
    },
    title='Percent Increase of Homers during the Full Moon, 2025',
)
fig.update_xaxes(tickformat='.1%')

import math
import os
import random
import re

import pandas as pd
import polars as pl
from tqdm.autonotebook import tqdm
from uncertainties import ufloat

from blob_cache import blob_cache


def rate_error(col, col_count):
    return 1.96 * (pl.col(col) * (1 - pl.col(col)) / pl.col(col_count)).sqrt()


def pipe_rate_error(_df: pl.LazyFrame | pl.DataFrame, col, col_count=None, out_name=None):
    if not out_name:
        out_name = f'{col}_error'
    if not col_count:
        col_count = f'{col}_count'
    if col_count not in _df.columns:
        _df = _df.with_columns(pl.col(col).count().alias(col_count))
    return _df.with_columns(rate_error(col, col_count).alias(out_name))


def make_fake_phase(_df: pl.LazyFrame, period: float = 26.7) -> pl.LazyFrame:

    new_data = (
        (
            pl
            .DataFrame(pl.date_range(pl.date(1900, 1, 1), pl.date(2026, 1, 1), eager=True))
            .lazy()
            .with_columns(days_since=(pl.col('date') - pl.col('date').min()).dt.total_days())
            .unique(subset=['date'])
            .select(['date', 'days_since'])
            .with_columns(
                phase=100 * (0.5 + 0.5 * (math.pi * pl.col('days_since') / period).cos())
            )
        )
        .sort('days_since')
        .with_columns(full_moon_flag=pl.col('phase').peak_max())
    )
    return _df.drop(['phase', 'full_moon_flag']).join(new_data, on='date', how='left')


def make_one_df_range(
    _merged,
    threshold=7,
    phase=99,
    range_size=1,
    fake_moon_period=None,
    group_by=['Moon Status'],
    stat='large_victory',
    max_full_moon=True,
):
    merged = _merged.lazy()
    if stat == 'large_victory':  # really should move stat creation out of this function
        merged = merged.with_columns(large_victory=pl.col('margin_of_victory').ge(threshold))
    if fake_moon_period:
        merged = merged.pipe(make_fake_phase, fake_moon_period)
    return (
        merged
        .with_columns(
            pl
            .when(max_full_moon)
            .then(pl.col('full_moon_flag'))
            .otherwise(pl.col('phase').is_between(phase, phase + range_size))
            .alias('Moon Status'),
        )
        .group_by(group_by, maintain_order=True)
        .agg(
            pl.col(stat).mean().alias(f'{stat}_rate'),
            pl.col('phase').mean(),
            pl.col('gid').count().alias('count'),
            pl.col(stat).sum().alias(f'{stat}_sum'),
            #    pl.col('daynight').eq('day').mean().alias('day_game_flag'),
        )
        .with_columns(blowout_threshold=threshold)
        .collect()
        .pipe(pipe_rate_error, f'{stat}_rate', 'count')
        .with_columns(phase_bin=phase, fake_moon_period=fake_moon_period)
    )


def clean_name(x):
    return x.replace('_', ' ').title()


def find_col(obj, pattern, regex=False):
    r"""Simple match column/key names in a dataframe/dictionary, optional regex. Without regex
    does as simple "s in col_name"

    usage:

    find_col(_df,'appraisal')

    or

    find_col(_df,'\w+_\w+',regex=True)

    """
    pattern = str(pattern)
    if not hasattr(obj, 'columns'):
        return []
    if not regex:
        return sorted([col for col in obj.columns if pattern in str(col)])
    else:
        return sorted([col for col in obj.columns if re.match(pattern, str(col))])


def make_phase_bin_data(
    merged,
    threshold=7,
    fake_moon_period=None,
    phase_range=[99],
    range_size=1,
    stat='large_victory',
    max_full_moon=True,
):

    phase_vs_all = (
        (
            pl.concat([
                make_one_df_range(
                    merged,
                    threshold=threshold,
                    phase=phase,
                    range_size=range_size,
                    fake_moon_period=fake_moon_period,
                    stat=stat,
                    max_full_moon=max_full_moon,
                )
                for phase in phase_range
            ]).pivot(on='Moon Status', index='phase_bin', values=[f'{stat}_rate', 'count'])
        )
        .pipe(pipe_rate_error, f'{stat}_rate_true', 'count_true')
        .pipe(pipe_rate_error, f'{stat}_rate_false', 'count_false')
    )
    phase_data_pandas = phase_vs_all.to_pandas()
    phase_data_pandas[f'{stat}_rate_in_phase_range'] = phase_data_pandas[
        [f'{stat}_rate_true', f'{stat}_rate_true_error']
    ].apply(lambda x: ufloat(x.iloc[0], x.iloc[1]), axis=1)
    phase_data_pandas[f'{stat}_rate_out_of_phase_range'] = phase_data_pandas[
        [f'{stat}_rate_false', f'{stat}_rate_false_error']
    ].apply(lambda x: ufloat(x.iloc[0], x.iloc[1]), axis=1)
    phase_data_pandas = phase_data_pandas.assign(
        ratio_ufloat=pd.col(f'{stat}_rate_in_phase_range')
        / pd.col(f'{stat}_rate_out_of_phase_range')
    ).assign(
        ratio=pd.col('ratio_ufloat').apply(lambda x: x.n),
        ratio_error=pd.col('ratio_ufloat').apply(lambda x: x.s),
        fake_moon_period=fake_moon_period,
    )
    return phase_data_pandas


def get_all_baseballdata(start_year=1900):
    return (
        pl
        .scan_csv(
            '/Volumes/BMD/Large Data Sets/baseball_data/gameinfo.csv',
            infer_schema_length=100000,
        )
        .with_columns(pl.col('date').cast(pl.String).str.to_date('%Y%m%d'))
        .filter(pl.col('date').gt(pl.date(start_year, 1, 1)))
        .with_columns(margin_of_victory=(pl.col('vruns') - pl.col('hruns')).abs())
    ).collect()


def make_moon(file='horizons_results.txt'):
    pandas_df = pd.read_table(file, sep=r'\s{4,}', parse_dates=True, engine='python')
    pandas_df.columns = ['date', 'phase']
    pandas_df['date'] = pd.to_datetime(pandas_df['date'])
    return (
        pl
        .from_pandas(pandas_df)
        .with_columns(pl.col('date').dt.date())
        .with_columns(full_moon_flag=pl.col('phase').peak_max())
    )


def get_merged_data():
    all_data = get_all_baseballdata()
    moon_data = make_moon()
    return all_data.join(moon_data, how='left', on=['date'])


@blob_cache(key='fake-moon-data2', keys_to_hash=['merged'])
def make_fake_vs_real_full_moons(
    merged,
    n=100,
    threshold=7,
    group=['Moon Status'],
    stat='large_victory',
    max_full_moon=True,
):
    out = []

    for _ in tqdm(range(n)):
        random_period = random.randint(250000, 350000) / 10000
        out.append(
            make_one_df_range(
                merged,
                fake_moon_period=random_period,
                threshold=threshold,
                group_by=group,
                stat=stat,
                max_full_moon=max_full_moon,
            )
        )
    out.append(
        make_one_df_range(
            merged,
            group_by=group,
            stat=stat,
            fake_moon_period=None,
            threshold=7,
        )
    )
    return (
        pl
        .concat(out)
        .pivot(
            on='Moon Status',
            index=['phase_bin', 'fake_moon_period']
            + (['season'] if 'season' in group else []),
            values=[f'{stat}_rate'],
            aggregate_function='first',
        )
        .with_columns(full_moon_ratio=pl.col('true') / 'false')
        .with_columns(pl.col('fake_moon_period').is_not_null().alias('fake moon'))
    )


def rename_for_gt(label):
    if label == 'Moon Status':
        return label
    label_kind = label.split('_')[-1]
    stat_name = '_'.join(label.split('_')[:-1])
    match label_kind:
        case 'count':
            return 'Game Count'
        case 'rate':
            return clean_name(label)
        case 'sum':
            return f'{clean_name(stat_name)} Game Count'
        case 'error':
            return clean_name(label)


def make_table(df, stat_name):
    return (
        df
        .drop(['fake_moon_period', 'blowout_threshold', 'phase', 'phase_bin'], strict=False)
        .sort('Moon Status')
        .style.fmt_percent(columns=[f'{stat_name}_rate', f'{stat_name}_rate_error'])
        .cols_label_with(fn=rename_for_gt)
        .cols_move(f'{stat_name}_rate_error', after=f'{stat_name}_rate')
        .tab_header(title=f'{clean_name(stat_name)} Full Moon Statistics')
    )