Bimodal Work Hours Distribution

import matplotlib.pyplot as plt
import numpy as np
from scipy import stats

np.random.seed(48)

# Bimodal: part-time and full-time workers
part_time = np.random.normal(4, 1, 300)
full_time = np.random.normal(8.5, 1.2, 700)
hours = np.concatenate([part_time, full_time])
hours = hours[(hours > 0) & (hours < 14)]

kde = stats.gaussian_kde(hours)
x = np.linspace(0, 14, 500)
y = kde(x)

fig, ax = plt.subplots(figsize=(12, 6), facecolor='#0a0a0f')
ax.set_facecolor('#0a0a0f')

# Neon cyan gradient
ax.fill_between(x, y, alpha=0.4, color='#27D3F5')
ax.plot(x, y, color='#27D3F5', linewidth=3)

# Add glow effect
ax.plot(x, y, color='#27D3F5', linewidth=6, alpha=0.3)

# Highlight modes
peaks = [4, 8.5]
for peak in peaks:
    ax.axvline(peak, color='#F5B027', linestyle='--', linewidth=2, alpha=0.8)
    kde_val = kde(peak)[0]
    ax.scatter([peak], [kde_val], color='#F5B027', s=100, zorder=5)
    ax.annotate(f'{peak}h', (peak, kde_val), textcoords='offset points', 
                xytext=(10, 10), color='#F5B027', fontsize=11, fontweight='bold')

ax.set_xlabel('Work Hours per Day', fontsize=12, color='white', fontweight='500')
ax.set_ylabel('Density', fontsize=12, color='white', fontweight='500')
ax.set_title('Daily Work Hours Distribution (Part-time vs Full-time)', fontsize=16, color='white', fontweight='bold', pad=15)

ax.tick_params(colors='white', labelsize=10)
for spine in ax.spines.values():
    spine.set_color('#334155')
ax.grid(True, alpha=0.1, color='white')
ax.set_xlim(0, 14)

plt.tight_layout()
plt.show()