import numpy as np from scipy.stats import chi2 from collections import Counter import pandas as pd from itertools import product def chi2_randomness(sequence, s, alpha=0.05): """ Тест Хи-квадрат на наличие марковской памяти в дискретной последовательности Args: sequence: Дискретная последовательность s: Порядок марковской цепи (глубина памяти) alpha: Уровень значимости Returns: Dict: Статистика, критическое значение, p-value и т.д. """ n = len(sequence) # количество элементов в исследуемой последовательности if n <= s: return { "statistic": 0.0, "critical_value": 0.0, "p_value": 1.0, "is_rejected": False, "df": 0, "warning": "Sequence too short" } alphabet = sorted(list(set(sequence))) # автоматически определяет алфавит системы (множество E) N = len(alphabet) # количество элементов в E # 1. Считаем только частоты (s+1)-цепочек (v_ij), N^{s+1} комбинаций v_ij_counts = Counter() # создаем обьект Counter for t in range(n - s): # проходим по всей истории окном длиной s+1. Каждая найденная комбинация сохраняется как кортеж (tuple) chain_s_plus_1 = tuple(sequence[t : t + s + 1]) # создаем ключ из (s+1)-цепочек v_ij_counts[chain_s_plus_1] += 1 # Если паттерн встретился впервые обьект Counter создаст новую запись # 2. Чтобы получить частоту v_i суммируем частоты v_ij по всем j для данного префикса i v_i_from_v_ij = Counter() # создаем обьект Counter for chain_s_plus_1, count in v_ij_counts.items(): prefix_i = chain_s_plus_1[:-1] # создаем новый ключ , берем первые s элементов v_i_from_v_ij[prefix_i] += count # считаем частоту v_i # 3. Вычисляем статистику хи-квадрат chi_square_stat = 0 # Проходим в цикле по всем встреченным s-цепочкам (i) for i_chain, v_i in v_i_from_v_ij.items(): expected = v_i / N # Ожидаемое значение согласно H0: v_i / N (равновероятность) # Для каждой i-цепочки смотрим все возможные продолжения j из алфавита (множество E) for j_state in alphabet: ij_chain = i_chain + (j_state,) v_ij = v_ij_counts.get(ij_chain, 0) chi_square_stat += ((v_ij - expected)**2) / expected df = (N**s) * (N - 1) # Степени свободы: N^s(N - 1) critical_value = chi2.ppf(1 - alpha, df) p_value = 1 - chi2.cdf(chi_square_stat, df) return { "statistic": chi_square_stat, "critical_value": critical_value, "p_value": p_value, "is_rejected": chi_square_stat > critical_value, "df": df } def quantile_transform(data, n_states=4)-> np.ndarray: """ Квантильное преобразование разбивает вещественные данные на N дискретных состояний на основе эмпирических квантилей. Это нейтрализует влияние тренда и обеспечивает P{E(j)} = 1/N Args: data: Одномерный массив вещественных чисел (приращения, доходности) n_states: Количество дискретных состояний (размер алфавита) Returns: np.ndarray: Последовательность целых чисел (индексов состояний) от 0 до N-1 """ # Рассчитываем пороги для N состояний # Для n_states=4 пороги будут на 25%, 50% и 75% quantiles = np.linspace(0, 100, n_states + 1)[1:-1] thresholds = np.percentile(data, quantiles) # Превращаем данные в последовательность состояний 0, 1, ..., N-1 return np.digitize(data, thresholds) def analyze_patterns(sequence, s_order): """ Анализирует отклонение частот цепочек длиной s+1 от теоретического ожидания """ n = len(sequence) s_plus_1 = s_order + 1 if n < s_plus_1: print(f"Warning: Sequence length ({n}) is too short for order s={s_order}") return pd.DataFrame() alphabet = sorted(list(set(sequence))) N = len(alphabet) # 1. Собираем все цепочки длины s+1 chains = [tuple(sequence[t : t + s_plus_1]) for t in range(n - s_order)] # список кортежей counts = Counter(chains) # Counter превращает список в словарь, где ключ — паттерн, а значение — сколько раз он встретился total_observed = len(chains) # 2. Ожидание для идеального хаоса num_combinations = N**s_plus_1 expected_count = total_observed / num_combinations # Evi_j = n/N^(s+1) results = [] # 3. Генерируем декартово произведение алфавита (все возможные цепочки состояний) for pattern in product(alphabet, repeat=s_plus_1): # например E x E x E = E^3 actual = counts.get(pattern, 0) # Считаем отклонение в процентах от ожидаемого при H0 diff_pct = ((actual - expected_count) / expected_count) * 100 if expected_count > 0 else 0 results.append({ "Pattern": pattern, "Actual": actual, "Expected": round(expected_count, 1), "Deviation_%": round(diff_pct, 2) }) return pd.DataFrame(results) if __name__ == "__main__": # --- 1: приращения цены с положительным матожиданием (тренд) --- print("\n1. приращения цены с положительным матожиданием (тренд):") data = np.random.normal(loc=0.5, scale=1.0, size=1000) # ТЕСТ 1: Разбиение по нулевому порогу seq_raw = (data > 0).astype(int) res_raw = chi2_randomness(seq_raw, s=1) # ТЕСТ 2: Разбиение по медиане seq_quant = quantile_transform(data, n_states=2) res_quant = chi2_randomness(seq_quant, s=1) counts_q = Counter(seq_quant) print(f" Без квантилей : X2={res_raw['statistic']:.2f}, P-value={res_raw['p_value']:.4f}") print(f" С квантилями (баланс {counts_q[0]}/{counts_q[1]}): X2={res_quant['statistic']:.2f}, P-value={res_quant['p_value']:.4f}") # --- 2: Модель марковской цепи маятник) --- print("\n2. модель марковской цепи маятник :") markov_seq = [0] for _ in range(999): curr = markov_seq[-1] # Вероятность сменить состояние - 90% markov_seq.append(np.random.choice([0, 1], p=[0.1, 0.9]) if curr == 0 else np.random.choice([0, 1], p=[0.9, 0.1])) res_markov = chi2_randomness(markov_seq, s=1) counts = Counter(markov_seq) print(f" Распределение: 0: {counts[0]}, 1: {counts[1]}") print(f" Результат теста: X2={res_markov['statistic']:.2f}, P-value={res_markov['p_value']:.4f}")