<?php
/**
 * The churn model: a small logistic regression trained on the store's own history.
 *
 * Training looks back one horizon: it takes every customer as they were on a past cutoff
 * date, and labels whether they bought again in the horizon after it. The model learns
 * which recency, frequency, spend and rhythm patterns came before customers leaving.
 * A held-out fifth of customers measures how well it works, against the simple rule, and
 * the model is only used when it beats the rule. Nothing leaves the site.
 *
 * @package Statbird
 */

namespace Statbird;

defined( 'ABSPATH' ) || exit;

/**
 * Feature building, training and prediction.
 */
class Model {

	const HORIZON_DAYS = 90;
	const MIN_SAMPLES  = 300;
	const MIN_CLASS    = 30;

	/**
	 * Feature vector for one customer as of a date.
	 *
	 * @param float $recency  Days since last purchase.
	 * @param int   $orders   Purchases so far.
	 * @param float $spent    Net spend so far.
	 * @param float $gap      Average days between purchases (0 for one purchase).
	 * @param float $tenure   Days since first purchase.
	 * @param float $expected Expected days between purchases (own gap, or the store's typical gap).
	 * @return float[]
	 */
	public static function features( $recency, $orders, $spent, $gap, $tenure, $expected ) {
		$recency  = max( 0.0, (float) $recency );
		$expected = max( 1.0, (float) $expected );
		return array(
			log1p( $recency ),
			log1p( max( 0, (int) $orders - 1 ) ),
			log1p( max( 0.0, (float) $spent ) ),
			log1p( max( 0.0, (float) $gap ) ),
			log1p( max( 0.0, (float) $tenure ) ),
			log1p( $recency / $expected ),
			$orders > 1 ? 1.0 : 0.0,
		);
	}

	/**
	 * The rule used when there is not enough history for a model: how overdue a customer is
	 * compared with their usual rhythm. 1.5 times their usual gap is about a 50% risk.
	 *
	 * @param float $recency  Days since last purchase.
	 * @param float $expected Expected days between purchases.
	 * @return float Probability 0 to 1.
	 */
	public static function rule_probability( $recency, $expected ) {
		$ratio = max( 0.0, (float) $recency ) / max( 1.0, (float) $expected );
		return 1 / ( 1 + exp( -2.2 * ( $ratio - 1.5 ) ) );
	}

	/**
	 * Logistic function.
	 *
	 * @param float $z Linear score.
	 * @return float
	 */
	public static function sigmoid( $z ) {
		if ( $z < -35 ) {
			return 0.0;
		}
		if ( $z > 35 ) {
			return 1.0;
		}
		return 1 / ( 1 + exp( -$z ) );
	}

	/**
	 * Trains on samples and returns a model array.
	 *
	 * @param array[] $x      Feature vectors.
	 * @param int[]   $y      Labels, 1 = did not come back.
	 * @param int     $epochs Gradient descent passes.
	 * @param float   $rate   Learning rate.
	 * @param float   $l2     Ridge penalty.
	 * @return array{w:float[],b:float,mean:float[],std:float[]}
	 */
	public static function train( array $x, array $y, $epochs = 400, $rate = 0.3, $l2 = 0.001 ) {
		$n = count( $x );
		$d = count( $x[0] );
		// Standardise so every feature learns at the same pace.
		$mean = array_fill( 0, $d, 0.0 );
		$std  = array_fill( 0, $d, 0.0 );
		foreach ( $x as $row ) {
			for ( $j = 0; $j < $d; $j++ ) {
				$mean[ $j ] += $row[ $j ] / $n;
			}
		}
		foreach ( $x as $row ) {
			for ( $j = 0; $j < $d; $j++ ) {
				$std[ $j ] += ( $row[ $j ] - $mean[ $j ] ) ** 2 / $n;
			}
		}
		for ( $j = 0; $j < $d; $j++ ) {
			$sd        = sqrt( $std[ $j ] );
			$std[ $j ] = $sd > 0 ? $sd : 1.0;
		}
		$z = array();
		foreach ( $x as $i => $row ) {
			for ( $j = 0; $j < $d; $j++ ) {
				$z[ $i ][ $j ] = ( $row[ $j ] - $mean[ $j ] ) / $std[ $j ];
			}
		}

		$w = array_fill( 0, $d, 0.0 );
		$b = 0.0;
		for ( $e = 0; $e < $epochs; $e++ ) {
			$gw = array_fill( 0, $d, 0.0 );
			$gb = 0.0;
			for ( $i = 0; $i < $n; $i++ ) {
				$s = $b;
				for ( $j = 0; $j < $d; $j++ ) {
					$s += $w[ $j ] * $z[ $i ][ $j ];
				}
				$err = self::sigmoid( $s ) - $y[ $i ];
				for ( $j = 0; $j < $d; $j++ ) {
					$gw[ $j ] += $err * $z[ $i ][ $j ];
				}
				$gb += $err;
			}
			for ( $j = 0; $j < $d; $j++ ) {
				$w[ $j ] -= $rate * ( $gw[ $j ] / $n + $l2 * $w[ $j ] );
			}
			$b -= $rate * $gb / $n;
		}
		return array(
			'w'    => $w,
			'b'    => $b,
			'mean' => $mean,
			'std'  => $std,
		);
	}

	/**
	 * Probability that a customer does not come back within the horizon.
	 *
	 * @param array   $model From train().
	 * @param float[] $f     Feature vector.
	 * @return float
	 */
	public static function predict( array $model, array $f ) {
		$s = $model['b'];
		foreach ( $model['w'] as $j => $wj ) {
			$s += $wj * ( $f[ $j ] - $model['mean'][ $j ] ) / $model['std'][ $j ];
		}
		return self::sigmoid( $s );
	}

	/**
	 * Area under the ROC curve: the chance a customer who left is scored riskier than one
	 * who stayed. 0.5 is guessing, 1 is perfect.
	 *
	 * @param float[] $scores Predicted probabilities.
	 * @param int[]   $labels 1 = left.
	 * @return float
	 */
	public static function auc( array $scores, array $labels ) {
		$pairs = array();
		foreach ( $scores as $i => $s ) {
			$pairs[] = array( $s, $labels[ $i ] );
		}
		usort(
			$pairs,
			static function ( $a, $b ) {
				return $a[0] <=> $b[0];
			}
		);
		$pos = array_sum( $labels );
		$neg = count( $labels ) - $pos;
		if ( ! $pos || ! $neg ) {
			return 0.5;
		}
		// Rank-sum with ties sharing their average rank.
		$rank_sum = 0.0;
		$n        = count( $pairs );
		for ( $i = 0; $i < $n; ) {
			$j = $i;
			while ( $j + 1 < $n && $pairs[ $j + 1 ][0] === $pairs[ $i ][0] ) {
				++$j;
			}
			$avg = ( $i + $j + 2 ) / 2;
			for ( $k = $i; $k <= $j; $k++ ) {
				if ( 1 === $pairs[ $k ][1] ) {
					$rank_sum += $avg;
				}
			}
			$i = $j + 1;
		}
		return ( $rank_sum - $pos * ( $pos + 1 ) / 2 ) / ( $pos * $neg );
	}

	/**
	 * Builds the training set from order history and trains, evaluating on a held-out fifth.
	 *
	 * @param float $store_gap The store's typical days between purchases.
	 * @param int   $now       Unix time to treat as today.
	 * @return array Model summary. 'use' says whether the model beat the rule.
	 */
	public static function fit( $store_gap, $now ) {
		global $wpdb;
		$h      = self::HORIZON_DAYS * DAY_IN_SECONDS;
		$cutoff = gmdate( 'Y-m-d H:i:s', $now - $h );
		$sql    = 'SELECT ckey,
				SUM(net > 0 AND created < %s) AS orders,
				SUM(CASE WHEN created < %s THEN GREATEST(net, 0) ELSE 0 END) AS spent,
				MIN(CASE WHEN net > 0 AND created < %s THEN created END) AS first_order,
				MAX(CASE WHEN net > 0 AND created < %s THEN created END) AS last_order,
				MAX(net > 0 AND created >= %s) AS came_back
			FROM ( ' . Orders::order_rows_sql() . ' ) t
			GROUP BY ckey
			HAVING orders > 0';
		// phpcs:ignore WordPress.DB.DirectDatabaseQuery.DirectQuery, WordPress.DB.DirectDatabaseQuery.NoCaching, WordPress.DB.PreparedSQL.NotPrepared -- inner SQL is built from constants.
		$rows = (array) $wpdb->get_results( $wpdb->prepare( $sql, $cutoff, $cutoff, $cutoff, $cutoff, $cutoff ), ARRAY_A );

		$x     = array();
		$y     = array();
		$rules = array();
		$at    = $now - $h;
		foreach ( $rows as $r ) {
			$first    = strtotime( $r['first_order'] . ' UTC' );
			$last     = strtotime( $r['last_order'] . ' UTC' );
			$orders   = (int) $r['orders'];
			$recency  = ( $at - $last ) / DAY_IN_SECONDS;
			$tenure   = ( $at - $first ) / DAY_IN_SECONDS;
			$gap      = $orders > 1 ? ( $last - $first ) / DAY_IN_SECONDS / ( $orders - 1 ) : 0;
			$expected = $orders > 1 ? max( 1, $gap ) : $store_gap;
			$x[]      = self::features( $recency, $orders, (float) $r['spent'], $gap, $tenure, $expected );
			$y[]      = (int) $r['came_back'] ? 0 : 1;
			$rules[]  = self::rule_probability( $recency, $expected );
		}

		// Large stores: a deterministic sample of 20,000 customers trains just as well, much faster.
		if ( count( $x ) > 20000 ) {
			$step = count( $x ) / 20000;
			$keep = array();
			for ( $k = 0; $k < 20000; $k++ ) {
				$keep[] = (int) floor( $k * $step );
			}
			$x     = array_values( array_intersect_key( $x, array_flip( $keep ) ) );
			$y     = array_values( array_intersect_key( $y, array_flip( $keep ) ) );
			$rules = array_values( array_intersect_key( $rules, array_flip( $keep ) ) );
		}
		$n   = count( $x );
		$pos = array_sum( $y );
		$out = array(
			'use'      => false,
			'samples'  => $n,
			'left'     => $pos,
			'horizon'  => self::HORIZON_DAYS,
			'trained'  => gmdate( 'Y-m-d H:i:s', $now ),
			'auc'      => null,
			'rule_auc' => null,
			'recall'   => null,
			'model'    => null,
		);
		if ( $n < self::MIN_SAMPLES || $pos < self::MIN_CLASS || ( $n - $pos ) < self::MIN_CLASS ) {
			$out['reason'] = 'not_enough_history';
			return $out;
		}

		// Deterministic 80/20 split by position, so results are repeatable.
		$tx = array();
		$ty = array();
		$hx = array();
		$hy = array();
		$hr = array();
		foreach ( $x as $i => $f ) {
			if ( 4 === $i % 5 ) {
				$hx[] = $f;
				$hy[] = $y[ $i ];
				$hr[] = $rules[ $i ];
			} else {
				$tx[] = $f;
				$ty[] = $y[ $i ];
			}
		}
		$model  = self::train( $tx, $ty );
		$scores = array();
		foreach ( $hx as $f ) {
			$scores[] = self::predict( $model, $f );
		}
		$auc      = self::auc( $scores, $hy );
		$rule_auc = self::auc( $hr, $hy );
		// Of the held-out customers who did leave, how many the model flagged at 50% or more.
		$flagged = 0;
		$left    = 0;
		foreach ( $hy as $i => $label ) {
			if ( 1 === $label ) {
				++$left;
				$flagged += $scores[ $i ] >= 0.5 ? 1 : 0;
			}
		}
		$out['auc']      = round( $auc, 3 );
		$out['rule_auc'] = round( $rule_auc, 3 );
		$out['recall']   = $left ? round( $flagged / $left, 3 ) : null;
		$out['model']    = $model;
		$out['use']      = $auc >= $rule_auc && $auc >= 0.6;
		if ( ! $out['use'] ) {
			$out['reason'] = 'rule_was_better';
		}
		// Retrain on everything for live predictions once the evaluation is recorded.
		if ( $out['use'] ) {
			$out['model'] = self::train( $x, $y );
		}
		return $out;
	}
}