Extended Poia Theory Simulation Framework
Extended Poia Theory Simulation Framework
The simulation allows you to:
- Test the core hypotheses of the Poia Theory across multiple domains
- Explore parameter sensitivities to optimize experimental designs
- Generate synthetic datasets that can guide statistical analysis approaches
- Visualize complex relationships between consciousness, quantum effects, and field interactions
"""
Comprehensive Poia Theory Simulation Framework
This program simulates the key experiments proposed to test the Poia Theory:
1. Neural-Quantum Correlation Experiment
2. Inter-Subject Neural Synchronization Experiment
3. Consciousness State Transitions and Field Effects
The simulation implements detailed models of neural dynamics, quantum systems,
and field interactions based on the mathematical formulations of the Poia Theory.
"""
import numpy as np
import matplotlib.pyplot as plt
from scipy import signal, stats
import pandas as pd
from sklearn.metrics import mutual_info_score
from tqdm import tqdm
import seaborn as sns
from mpl_toolkits.mplot3d import Axes3D
from matplotlib.animation import FuncAnimation
import networkx as nx
from scipy.integrate import solve_ivp
import warnings
warnings.filterwarnings('ignore')
class PoiaSimulation:
def __init__(self, params=None):
"""Initialize the comprehensive simulation framework."""
# Default parameters
self.params = {
# General parameters
'random_seed': 42, # Random seed for reproducibility
'verbose': True, # Print detailed information
# Neural parameters
'n_channels': 64, # Number of EEG channels
'fs': 500, # Sampling frequency (Hz)
'trial_duration': 15, # Trial duration in seconds
'n_trials': 200, # Number of trials per condition
'n_subjects': 24, # Number of subjects
# Brain regions and connectivity
'n_regions': 10, # Number of brain regions
'region_names': ['PFC', 'ACC', 'DLPFC', 'OFC', 'VC', 'TC', 'PC', 'HC', 'BG', 'TH'],
'connectivity_density': 0.3, # Density of connectivity matrix
# Coherence parameters
'phi_critical': 0.27, # Critical integrated information threshold
'phi_std': 0.04, # Standard deviation of phi across subjects
'alpha_band': [8, 12], # Alpha frequency band (Hz)
'theta_band': [4, 7], # Theta frequency band (Hz)
'gamma_band': [30, 45], # Gamma frequency band (Hz)
# Quantum parameters
'lambda_0': 1e-6, # Base coupling parameter
'alpha': 1.5, # Scaling exponent for lambda(Phi)
'qrng_bias': 0.5, # Baseline QRNG bias (0.5 = unbiased)
'amplification': 1e3, # Amplification factor for quantum effects
'quantum_coherence_time': 1e-6, # Quantum coherence time (s)
# Field coupling parameters
'k_alpha': 0.18, # Alpha band coupling constant
'k_theta': 0.12, # Theta band coupling constant
'k_gamma': 0.05, # Gamma band coupling constant
'distance_decay': 2.0, # Exponent for distance decay (inverse square = 2.0)
'field_propagation_speed': 3e8, # Field propagation speed (m/s) - speed of light
# Inter-subject parameters
'distance_levels': [0.5, 2.0, 5.0, 10.0, 20.0], # Distance between subjects (m)
'baseline_psi': 0.05, # Baseline phase synchronization index
# Experimental conditions
'conditions': ['rest', 'focused', 'open', 'cognitive', 'control'],
'condition_phi_mod': { # Phi modification by condition
'rest': 1.0, # Baseline
'focused': 1.2, # Increased in focused attention
'open': 1.1, # Slightly increased in open monitoring
'cognitive': 0.9, # Decreased during cognitive load
'control': 0.0, # No effect in control condition
},
# Noise parameters
'neural_noise': 0.2, # Neural measurement noise level
'quantum_noise': 0.01, # Quantum measurement noise level
'em_noise': 0.05, # Environmental EM noise level
# Syntropy parameters
'eta_critical': 0.2, # Critical efficiency parameter
'system_size_scaling': 0.5, # Scaling exponent for system size (alpha)
'coupling_scaling': 0.8, # Scaling exponent for coupling strength (beta)
}
# Update with custom parameters if provided
if params:
self.params.update(params)
# Set random seed
np.random.seed(self.params['random_seed'])
# Initialize data structures
self.subjects = self._initialize_subjects()
self.brain_connectivity = self._initialize_brain_connectivity()
self.results = {
'neural_quantum': {},
'inter_subject': {},
'consciousness_states': {},
'syntropy': {}
}
def _initialize_subjects(self):
"""Initialize subject parameters with individual variations."""
subjects = []
for i in range(self.params['n_subjects']):
# Generate subject-specific parameters with natural variation
phi_base = np.random.normal(
self.params['phi_critical'],
self.params['phi_std']
)
# Ensure phi_base is positive
phi_base = max(0.1, phi_base)
# Add individual frequency peaks with variation
alpha_peak = np.random.normal(10, 0.5) # Center of alpha band
theta_peak = np.random.normal(5.5, 0.4) # Center of theta band
gamma_peak = np.random.normal(38, 2.0) # Center of gamma band
# Regional permittivity and conductivity variations
region_permittivity = {}
region_conductivity = {}
for region in self.params['region_names']:
# Base values with individual variation
if region in ['PFC', 'DLPFC', 'OFC']:
# Prefrontal regions
perm_factor = np.random.normal(1.15, 0.12)
cond_factor = np.random.normal(0.92, 0.11)
elif region in ['VC']:
# Visual cortex
perm_factor = np.random.normal(0.88, 0.09)
cond_factor = np.random.normal(1.06, 0.13)
elif region in ['HC']:
# Hippocampus
perm_factor = np.random.normal(1.32, 0.15)
cond_factor = np.random.normal(0.95, 0.12)
else:
# Other regions
perm_factor = np.random.normal(1.0, 0.1)
cond_factor = np.random.normal(1.0, 0.1)
region_permittivity[region] = perm_factor
region_conductivity[region] = cond_factor
subject = {
'id': i,
'phi_base': phi_base,
'alpha_peak': alpha_peak,
'theta_peak': theta_peak,
'gamma_peak': gamma_peak,
'neural_noise_factor': np.random.uniform(0.8, 1.2),
'quantum_susceptibility': np.random.uniform(0.7, 1.3),
'region_permittivity': region_permittivity,
'region_conductivity': region_conductivity,
'consciousness_state': 'awake', # Initial state
'eta': np.random.normal(0.6, 0.1) # Syntropy efficiency parameter
}
subjects.append(subject)
return subjects
def _initialize_brain_connectivity(self):
"""Initialize brain connectivity matrix between regions."""
n_regions = self.params['n_regions']
density = self.params['connectivity_density']
# Create a random connectivity matrix with specified density
connectivity = np.random.rand(n_regions, n_regions)
threshold = np.percentile(connectivity, 100 * (1 - density))
connectivity[connectivity < threshold] = 0
# Make it symmetric (undirected graph)
connectivity = np.maximum(connectivity, connectivity.T)
# Set diagonal to zero (no self-connections)
np.fill_diagonal(connectivity, 0)
# Normalize
if np.max(connectivity) > 0:
connectivity = connectivity / np.max(connectivity)
return connectivity
def _generate_neural_oscillations(self, duration, condition, subject):
"""Generate simulated EEG data with realistic oscillations."""
fs = self.params['fs']
n_samples = int(duration * fs)
n_channels = self.params['n_channels']
n_regions = self.params['n_regions']
# Time vector
t = np.arange(n_samples) / fs
# Base oscillations (pink noise as background)
eeg_data = np.zeros((n_channels, n_samples))
for ch in range(n_channels):
# Generate pink noise (1/f)
pink_noise = self._generate_pink_noise(n_samples)
eeg_data[ch, :] = pink_noise
# Add oscillatory components with condition-specific modulation
condition_mod = self.params['condition_phi_mod'][condition]
# Skip oscillations for control condition
if condition != 'control':
# Alpha oscillations
alpha_amp = 1.0 * condition_mod
if condition == 'focused':
alpha_amp *= 1.3 # Enhanced alpha in focused attention
alpha_freq = subject['alpha_peak']
alpha_phase = np.random.uniform(0, 2*np.pi)
alpha_osc = alpha_amp * np.sin(2 * np.pi * alpha_freq * t + alpha_phase)
# Theta oscillations
theta_amp = 0.7 * condition_mod
if condition == 'open':
theta_amp *= 1.4 # Enhanced theta in open monitoring
theta_freq = subject['theta_peak']
theta_phase = np.random.uniform(0, 2*np.pi)
theta_osc = theta_amp * np.sin(2 * np.pi * theta_freq * t + theta_phase)
# Gamma oscillations
gamma_amp = 0.4 * condition_mod
if condition == 'cognitive':
gamma_amp *= 1.5 # Enhanced gamma in cognitive tasks
gamma_freq = subject['gamma_peak']
gamma_phase = np.random.uniform(0, 2*np.pi)
gamma_osc = gamma_amp * np.sin(2 * np.pi * gamma_freq * t + gamma_phase)
# Generate region-specific activity
region_activity = np.zeros((n_regions, n_samples))
for r in range(n_regions):
# Base oscillations for this region
region_name = self.params['region_names'][r]
# Region-specific frequency modulation
if region_name in ['PFC', 'ACC', 'DLPFC']:
# Frontal regions - enhanced theta
r_theta = theta_osc * 1.2
r_alpha = alpha_osc * 0.8
r_gamma = gamma_osc * 0.9
elif region_name in ['VC', 'TC']:
# Sensory regions - enhanced alpha
r_theta = theta_osc * 0.7
r_alpha = alpha_osc * 1.3
r_gamma = gamma_osc * 1.1
elif region_name in ['HC']:
# Hippocampus - enhanced theta
r_theta = theta_osc * 1.4
r_alpha = alpha_osc * 0.6
r_gamma = gamma_osc * 0.8
else:
# Other regions - balanced
r_theta = theta_osc
r_alpha = alpha_osc
r_gamma = gamma_osc
# Combine oscillations
region_activity[r, :] = r_theta + r_alpha + r_gamma
# Add region-specific noise
noise_level = 0.1 * self.params['neural_noise']
region_activity[r, :] += noise_level * np.random.randn(n_samples)
# Propagate activity between regions based on connectivity
propagated_activity = np.zeros_like(region_activity)
for r1 in range(n_regions):
for r2 in range(n_regions):
if self.brain_connectivity[r1, r2] > 0:
# Add weighted contribution from connected regions
# With random delay based on connection strength
delay_samples = int(fs * 0.01 * (1 - self.brain_connectivity[r1, r2]))
if delay_samples > 0 and delay_samples < n_samples:
propagated_activity[r1, :-delay_samples] += (
self.brain_connectivity[r1, r2] *
region_activity[r2, delay_samples:]
)
# Combine original and propagated activity
region_activity = 0.6 * region_activity + 0.4 * propagated_activity
# Map region activity to EEG channels
# Each channel gets weighted contribution from regions
for ch in range(n_channels):
for r in range(n_regions):
# Simple spatial mapping - different channels get different weights from regions
# In a real implementation, this would use a head model and lead field matrix
weight = np.exp(-(ch/n_channels - r/n_regions)**2 / 0.1)
eeg_data[ch, :] += 0.5 * weight * region_activity[r, :]
# Add measurement noise
noise_level = self.params['neural_noise'] * subject['neural_noise_factor']
noise = noise_level * np.random.randn(n_channels, n_samples)
eeg_data += noise
return eeg_data, region_activity
def _generate_pink_noise(self, n_samples):
"""Generate pink noise (1/f spectrum)."""
white_noise = np.random.randn(n_samples)
# Create pink noise by filtering white noise
b, a = signal.butter(1, 0.9, 'low')
pink_noise = signal.lfilter(b, a, white_noise)
# Normalize
pink_noise = pink_noise / np.std(pink_noise)
return pink_noise
def _calculate_neural_coherence(self, eeg_data, band):
"""Calculate neural coherence in specified frequency band."""
fs = self.params['fs']
n_channels = self.params['n_channels']
# Apply bandpass filter
low, high = band
b, a = signal.butter(3, [low/(fs/2), high/(fs/2)], 'bandpass')
filtered_data = np.zeros_like(eeg_data)
for ch in range(n_channels):
filtered_data[ch, :] = signal.filtfilt(b, a, eeg_data[ch, :])
# Calculate Hilbert transform to get analytic signal
analytic_signal = signal.hilbert(filtered_data)
phase = np.angle(analytic_signal)
amplitude = np.abs(analytic_signal)
# Calculate phase synchronization index (PSI) across channels
n_pairs = 0
psi_sum = 0
for i in range(n_channels):
for j in range(i+1, n_channels):
phase_diff = phase[i, :] - phase[j, :]
psi = np.abs(np.mean(np.exp(1j * phase_diff)))
psi_sum += psi
n_pairs += 1
# Average PSI across all channel pairs
mean_psi = psi_sum / n_pairs if n_pairs > 0 else 0
# Also calculate power in this band
power = np.mean(amplitude**2)
return mean_psi, power
def _calculate_phi(self, eeg_data, region_activity, subject, condition):
"""Calculate integrated information (Phi) from neural data."""
# Calculate coherence in different frequency bands
alpha_coherence, alpha_power = self._calculate_neural_coherence(eeg_data, self.params['alpha_band'])
theta_coherence, theta_power = self._calculate_neural_coherence(eeg_data, self.params['theta_band'])
gamma_coherence, gamma_power = self._calculate_neural_coherence(eeg_data, self.params['gamma_band'])
# Calculate effective connectivity between regions
n_regions = self.params['n_regions']
effective_conn = np.zeros((n_regions, n_regions))
for i in range(n_regions):
for j in range(n_regions):
if i != j:
# Calculate Granger-like causality (simplified)
# In a real implementation, this would use proper Granger causality or transfer entropy
corr = np.corrcoef(region_activity[i, 1:], region_activity[j, :-1])[0, 1]
effective_conn[i, j] = max(0, corr)
# Calculate integration (I) - simplified measure based on connectivity
integration = np.mean(effective_conn)
# Calculate differentiation (D) - based on spectral diversity
band_powers = [alpha_power, theta_power, gamma_power]
differentiation = np.std(band_powers) / np.mean(band_powers)
# Calculate Phi as a combination of integration and differentiation
phi_raw = integration * differentiation
# Apply subject-specific baseline and condition modifier
condition_mod = self.params['condition_phi_mod'][condition]
phi = subject['phi_base'] * condition_mod * phi_raw
# Add some random variation
phi = max(0, phi + np.random.normal(0, 0.05))
# Additional metrics for comprehensive Phi calculation
metrics = {
'alpha_coherence': alpha_coherence,
'theta_coherence': theta_coherence,
'gamma_coherence': gamma_coherence,
'integration': integration,
'differentiation': differentiation,
'phi_raw': phi_raw,
'phi': phi
}
return phi, metrics
def _calculate_lambda(self, phi):
"""Calculate lambda(Phi) scaling function."""
phi_critical = self.params['phi_critical']
# If phi is below critical threshold, lambda is zero
if phi < phi_critical:
return 0
# Otherwise, calculate lambda based on scaling function
lambda_0 = self.params['lambda_0']
alpha = self.params['alpha']
lambda_phi = lambda_0 * ((phi - phi_critical) ** alpha)
return lambda_phi
def _simulate_quantum_system(self, duration, phi, subject):
"""Simulate quantum system with potential influence from phi."""
fs = self.params['fs']
n_samples = int(duration * fs)
# Calculate lambda based on phi
lambda_phi = self._calculate_lambda(phi)
# Apply subject-specific quantum susceptibility
lambda_phi *= subject['quantum_susceptibility']
# Apply amplification factor
effect_size = lambda_phi * self.params['amplification']
# Base probability (0.5 for unbiased)
p = self.params['qrng_bias']
# Time vector
t = np.arange(n_samples) / fs
# Simulate quantum coherence time effects
coherence_time = self.params['quantum_coherence_time']
n_coherence_periods = int(np.ceil(duration / coherence_time))
# Generate quantum random numbers with potential phi influence
qrng_output = np.zeros(n_samples)
for i in range(n_coherence_periods):
start_idx = int(i * coherence_time * fs)
end_idx = min(int((i + 1) * coherence_time * fs), n_samples)
if end_idx <= start_idx:
continue
# Modified probability based on phi
p_mod = p + effect_size
# Ensure probability is valid
p_mod = max(0, min(1, p_mod))
# Generate binary outcome for this coherence period
outcome = np.random.binomial(1, p_mod)
qrng_output[start_idx:end_idx] = outcome
# Add measurement noise
noise = self.params['quantum_noise'] * np.random.randn(n_samples)
qrng_output_noisy = qrng_output + noise
# Calculate bias
bias = np.mean(qrng_output)
return qrng_output_noisy, bias
def _calculate_field_coupling(self, subject1, subject2, distance, band):
"""Calculate field coupling strength between two subjects at given distance."""
# Get band-specific coupling constant
if band == 'alpha':
k_band = self.params['k_alpha']
elif band == 'theta':
k_band = self.params['k_theta']
elif band == 'gamma':
k_band = self.params['k_gamma']
else:
k_band = 0.1 # Default
# Calculate distance decay
distance_decay = self.params['distance_decay']
distance_factor = distance ** (-distance_decay)
# Calculate coupling strength
coupling = k_band * distance_factor
# Adjust for subject-specific factors
# Average permittivity and conductivity across regions
perm1 = np.mean(list(subject1['region_permittivity'].values()))
perm2 = np.mean(list(subject2['region_permittivity'].values()))
cond1 = np.mean(list(subject1['region_conductivity'].values()))
cond2 = np.mean(list(subject2['region_conductivity'].values()))
# Coupling depends on permittivity and inversely on conductivity
subject_factor = np.sqrt((perm1 * perm2) / (cond1 * cond2))
coupling *= subject_factor
return coupling
def _simulate_inter_subject_synchronization(self, subject1_idx, subject2_idx, distance, condition):
"""Simulate neural synchronization between two subjects at specified distance."""
subject1 = self.subjects[subject1_idx]
subject2 = self.subjects[subject2_idx]
duration = self.params['trial_duration']
# Generate independent neural activity for each subject
eeg_data1, region_activity1 = self._generate_neural_oscillations(duration, condition, subject1)
eeg_data2, region_activity2 = self._generate_neural_oscillations(duration, condition, subject2)
# Calculate Phi for each subject
phi1, metrics1 = self._calculate_phi(eeg_data1, region_activity1, subject1, condition)
phi2, metrics2 = self._calculate_phi(eeg_data2, region_activity2, subject2, condition)
# Calculate field coupling for different frequency bands
alpha_coupling = self._calculate_field_coupling(subject1, subject2, distance, 'alpha')
theta_coupling = self._calculate_field_coupling(subject1, subject2, distance, 'theta')
gamma_coupling = self._calculate_field_coupling(subject1, subject2, distance, 'gamma')
# Calculate expected synchronization based on coupling and Phi similarity
phi_similarity = np.exp(-abs(phi1 - phi2) / 0.3) # Similarity decreases with Phi difference
# Baseline synchronization (chance level)
baseline_psi = self.params['baseline_psi']
# Calculate expected PSI for each band
alpha_psi_expected = baseline_psi + alpha_coupling * phi_similarity
theta_psi_expected = baseline_psi + theta_coupling * phi_similarity
gamma_psi_expected = baseline_psi + gamma_coupling * phi_similarity
# Add noise to synchronization
alpha_psi = max(0, min(1, alpha_psi_expected + np.random.normal(0, 0.02)))
theta_psi = max(0, min(1, theta_psi_expected + np.random.normal(0, 0.02)))
gamma_psi = max(0, min(1, gamma_psi_expected + np.random.normal(0, 0.02)))
# Calculate information transfer between subjects
# In a real implementation, this would use transfer entropy or similar measures
# Here we use a simplified model based on synchronization and Phi
transfer_efficiency = 0.05 * (alpha_psi + theta_psi + gamma_psi) * np.sqrt(phi1 * phi2)
# Add noise to transfer efficiency
transfer_efficiency = max(0, transfer_efficiency + np.random.normal(0, 0.01))
return {
'subject1': subject1_idx,
'subject2': subject2_idxdata, ax=ax2)
ax2.set_title('Phase Synchronization by Condition')
ax2.set_ylabel('PSI')
# 3. Plot transfer efficiency vs Phi similarity
ax3 = fig.add_subplot(2, 2, 3)
sns.scatterplot(x='phi_similarity', y='transfer_efficiency', hue='distance', data=df, ax=ax3)
ax3.set_title('Information Transfer vs Phi Similarity')
ax3.set_xlabel('Phi Similarity')
ax3.set_ylabel('Transfer Efficiency')
# 4. Plot distance scaling
ax4 = fig.add_subplot(2, 2, 4)
# Get scaling results
scaling = self.results['inter_subject']['distance_scaling']
# Plot theoretical inverse square law
x = np.logspace(np.log10(min(distances)), np.log10(max(distances)), 100)
y_inverse_square = 0.1 * x**(-2) + self.params['baseline_psi']
# Plot measured scaling
y_alpha = np.exp(scaling['alpha_intercept']) * x**(scaling['alpha_exponent']) + self.params['baseline_psi']
y_theta = np.exp(scaling['theta_intercept']) * x**(scaling['theta_exponent']) + self.params['baseline_psi']
y_gamma = np.exp(scaling['gamma_intercept']) * x**(scaling['gamma_exponent']) + self.params['baseline_psi']
ax4.plot(x, y_inverse_square, 'k--', label=f'Inverse Square Law (r^-2)')
ax4.plot(x, y_alpha, '-', label=f'Alpha (r^{scaling["alpha_exponent"]:.2f})')
ax4.plot(x, y_theta, '-', label=f'Theta (r^{scaling["theta_exponent"]:.2f})')
ax4.plot(x, y_gamma, '-', label=f'Gamma (r^{scaling["gamma_exponent"]:.2f})')
ax4.set_xscale('log')
ax4.set_yscale('log')
ax4.set_title('Distance Scaling of Synchronization')
ax4.set_xlabel('Distance (m)')
ax4.set_ylabel('PSI - Baseline')
ax4.legend()
plt.tight_layout()
return fig
def plot_consciousness_transitions(self):
"""Plot consciousness state transition results."""
if 'transitions' not in self.results['consciousness_states']:
print("No consciousness transition results to plot. Run the experiment first.")
return None
transitions = self.results['consciousness_states']['transitions']
# Select a few transitions to plot
if len(transitions) > 3:
plot_transitions = [transitions[0], transitions[3], transitions[5]]
else:
plot_transitions = transitions
# Set up the figure
fig = plt.figure(figsize=(15, 12))
for i, transition in enumerate(plot_transitions):
# Extract data
time = transition['time']
phi = transition['phi_trajectory']
p_conscious = transition['p_conscious']
alpha_power = transition['alpha_power_trajectory']
theta_power = transition['theta_power_trajectory']
gamma_power = transition['gamma_power_trajectory']
# Plot Phi and consciousness probability
ax1 = fig.add_subplot(3, 3, i*3 + 1)
ax1.plot(time, phi, 'b-', label='Phi')
ax1.set_ylabel('Phi')
ax1.set_xlabel('Time (s)')
ax1_twin = ax1.twinx()
ax1_twin.plot(time, p_conscious, 'r--', label='P(Conscious)')
ax1_twin.set_ylabel('P(Conscious)')
ax1_twin.set_ylim(0, 1)
# Add horizontal line at Phi_critical
ax1.axhline(y=self.params['phi_critical'], color='k', linestyle=':', alpha=0.5)
# Add title
ax1.set_title(f"Subject {transition['subject']}: {transition['initial_state']} → {transition['target_state']}")
# Combine legends
lines1, labels1 = ax1.get_legend_handles_labels()
lines2, labels2 = ax1_twin.get_legend_handles_labels()
ax1.legend(lines1 + lines2, labels1 + labels2, loc='upper right')
# Plot power in different frequency bands
ax2 = fig.add_subplot(3, 3, i*3 + 2)
ax2.plot(time, alpha_power, 'g-', label='Alpha')
ax2.plot(time, theta_power, 'm-', label='Theta')
ax2.plot(time, gamma_power, 'c-', label='Gamma')
ax2.set_ylabel('Relative Power')
ax2.set_xlabel('Time (s)')
ax2.legend()
ax2.set_title('Frequency Band Power')
# Plot 3D trajectory
ax3 = fig.add_subplot(3, 3, i*3 + 3, projection='3d')
ax3.plot(alpha_power, theta_power, gamma_power, 'k-')
ax3.scatter(alpha_power[0], theta_power[0], gamma_power[0], c='g', s=50, label='Start')
ax3.scatter(alpha_power[-1], theta_power[-1], gamma_power[-1], c='r', s=50, label='End')
# Add color gradient based on time
points = ax3.scatter(
alpha_power, theta_power, gamma_power,
c=time, cmap='viridis',
s=10, alpha=0.5
)
ax3.set_xlabel('Alpha Power')
ax3.set_ylabel('Theta Power')
ax3.set_zlabel('Gamma Power')
ax3.set_title('State Space Trajectory')
ax3.legend()
# Add colorbar
cbar = plt.colorbar(points, ax=ax3, pad=0.1)
cbar.set_label('Time (s)')
plt.tight_layout()
# Create a second figure for transition point analysis
if 'transition_points' in self.results['consciousness_states']:
tp_df = self.results['consciousness_states']['transition_points']
fig2 = plt.figure(figsize=(12, 6))
# Plot Phi at transition points
ax_tp = fig2.add_subplot(1, 1, 1)
sns.boxplot(x='transition', y='phi', data=tp_df, ax=ax_tp)
ax_tp.set_title('Phi at Consciousness Transition Points')
ax_tp.set_xlabel('Transition')
ax_tp.set_ylabel('Phi')
# Add horizontal line at mean transition Phi
if 'mean_transition_phi' in self.results['consciousness_states']:
mean_phi = self.results['consciousness_states']['mean_transition_phi']
ax_tp.axhline(y=mean_phi, color='r', linestyle='--', alpha=0.7)
# Add text annotation
ax_tp.text(
0.02, 0.95, f'Mean Transition Phi: {mean_phi:.3f}',
transform=ax_tp.transAxes, fontsize=12,
verticalalignment='top', bbox=dict(boxstyle='round', facecolor='white', alpha=0.7)
)
# Add horizontal line at Phi_critical
ax_tp.axhline(y=self.params['phi_critical'], color='k', linestyle=':', alpha=0.5)
ax_tp.text(
0.02, 0.85, f'Phi Critical: {self.params["phi_critical"]:.3f}',
transform=ax_tp.transAxes, fontsize=12,
verticalalignment='top', bbox=dict(boxstyle='round', facecolor='white', alpha=0.7)
)
plt.xticks(rotation=45, ha='right')
plt.tight_layout()
return fig, fig2
return fig
def plot_syntropy_results(self):
"""Plot syntropy analysis results."""
if 'raw_data' not in self.results['syntropy']:
print("No syntropy results to plot. Run the analysis first.")
return None
df = self.results['syntropy']['raw_data']
# Set up the figure
fig = plt.figure(figsize=(15, 12))
# 1. Plot syntropic vs non-syntropic systems
ax1 = fig.add_subplot(2, 2, 1)
# Create a 2D histogram of system size vs coupling strength, colored by syntropic behavior
syntropic_df = df[df['is_syntropic']]
non_syntropic_df = df[~df['is_syntropic']]
ax1.scatter(
non_syntropic_df['system_size'],
non_syntropic_df['coupling_strength'],
c='blue', alpha=0.5, label='Non-Syntropic'
)
ax1.scatter(
syntropic_df['system_size'],
syntropic_df['coupling_strength'],
c='red', alpha=0.5, label='Syntropic'
)
ax1.set_xscale('log')
ax1.set_yscale('log')
ax1.set_xlabel('System Size (N)')
ax1.set_ylabel('Coupling Strength (J)')
ax1.set_title('Syntropic vs Non-Syntropic Systems')
ax1.legend()
# 2. Plot critical efficiency threshold
ax2 = fig.add_subplot(2, 2, 2)
# Create a 3D scatter plot
scatter = ax2.scatter(
df['system_size'],
df['coupling_strength'],
df['eta_critical'],
c=df['is_syntropic'], cmap='coolwarm', alpha=0.7
)
ax2.set_xscale('log')
ax2.set_yscale('log')
ax2.set_xlabel('System Size (N)')
ax2.set_ylabel('Coupling Strength (J)')
ax2.set_title('Critical Efficiency Threshold')
# Add colorbar
cbar = plt.colorbar(scatter, ax=ax2)
cbar.set_label('Syntropic (1) vs Non-Syntropic (0)')
# 3. Plot scaling law verification
ax3 = fig.add_subplot(2, 2, 3)
# Get scaling law results
scaling = self.results['syntropy']['scaling_law']
# Create predicted vs actual plot
log_size = np.log(df['system_size'])
log_coupling = np.log(df['coupling_strength'])
log_eta_actual = np.log(df['eta_critical'])
# Predicted values from scaling law
log_eta_pred = scaling['intercept'] + scaling['size_exponent'] * log_size + scaling['coupling_exponent'] * log_coupling
ax3.scatter(log_eta_actual, log_eta_pred, alpha=0.5)
# Add diagonal line
min_val = min(log_eta_actual.min(), log_eta_pred.min())
max_val = max(log_eta_actual.max(), log_eta_pred.max())
ax3.plot([min_val, max_val], [min_val, max_val], 'k--')
ax3.set_xlabel('Actual log(η_critical)')
ax3.set_ylabel('Predicted log(η_critical)')
ax3.set_title('Scaling Law Verification')
# Add text with scaling exponents
ax3.text(
0.05, 0.95,
f"Size exponent: {scaling['size_exponent']:.3f} (expected: {scaling['expected_size_exponent']:.3f})\n"
f"Coupling exponent: {scaling['coupling_exponent']:.3f} (expected: {scaling['expected_coupling_exponent']:.3f})",
transform=ax3.transAxes, fontsize=10,
verticalalignment='top', bbox=dict(boxstyle='round', facecolor='white', alpha=0.7)
)
# 4. Plot efficiency vs energy input
ax4 = fig.