|
| 1 | +""" pyplots.ai |
| 2 | +manhattan-gwas: Manhattan Plot for GWAS |
| 3 | +Library: plotly 6.5.0 | Python 3.13.11 |
| 4 | +Quality: 92/100 | Created: 2025-12-31 |
| 5 | +""" |
| 6 | + |
| 7 | +import numpy as np |
| 8 | +import pandas as pd |
| 9 | +import plotly.graph_objects as go |
| 10 | + |
| 11 | + |
| 12 | +# Data - Simulated GWAS results |
| 13 | +np.random.seed(42) |
| 14 | + |
| 15 | +# Chromosome lengths (simplified, in Mb) |
| 16 | +chr_lengths = { |
| 17 | + "1": 249, |
| 18 | + "2": 243, |
| 19 | + "3": 198, |
| 20 | + "4": 191, |
| 21 | + "5": 182, |
| 22 | + "6": 171, |
| 23 | + "7": 159, |
| 24 | + "8": 146, |
| 25 | + "9": 141, |
| 26 | + "10": 136, |
| 27 | + "11": 135, |
| 28 | + "12": 134, |
| 29 | + "13": 115, |
| 30 | + "14": 107, |
| 31 | + "15": 102, |
| 32 | + "16": 90, |
| 33 | + "17": 83, |
| 34 | + "18": 80, |
| 35 | + "19": 59, |
| 36 | + "20": 64, |
| 37 | + "21": 47, |
| 38 | + "22": 51, |
| 39 | +} |
| 40 | + |
| 41 | +# Generate SNPs for each chromosome |
| 42 | +data = [] |
| 43 | +cumulative_pos = 0 |
| 44 | +chr_centers = {} |
| 45 | + |
| 46 | +for chrom, length in chr_lengths.items(): |
| 47 | + # Number of SNPs proportional to chromosome length |
| 48 | + n_snps = int(length * 40) |
| 49 | + positions = np.sort(np.random.uniform(0, length * 1e6, n_snps)) |
| 50 | + |
| 51 | + # Generate p-values (mostly non-significant, with some peaks) |
| 52 | + pvalues = np.random.uniform(0, 1, n_snps) |
| 53 | + |
| 54 | + # Add significant peaks on chromosomes 2, 8, and 15 |
| 55 | + if chrom == "2": |
| 56 | + peak_idx = np.abs(positions - 100e6).argmin() |
| 57 | + pvalues[peak_idx - 5 : peak_idx + 5] = 10 ** (-np.random.uniform(8, 12, 10)) |
| 58 | + elif chrom == "8": |
| 59 | + peak_idx = np.abs(positions - 70e6).argmin() |
| 60 | + pvalues[peak_idx - 3 : peak_idx + 3] = 10 ** (-np.random.uniform(7.5, 10, 6)) |
| 61 | + elif chrom == "15": |
| 62 | + peak_idx = np.abs(positions - 50e6).argmin() |
| 63 | + pvalues[peak_idx - 4 : peak_idx + 4] = 10 ** (-np.random.uniform(9, 14, 8)) |
| 64 | + |
| 65 | + # Calculate cumulative position |
| 66 | + cumulative_positions = positions + cumulative_pos |
| 67 | + chr_centers[chrom] = cumulative_pos + (length * 1e6) / 2 |
| 68 | + |
| 69 | + for i in range(n_snps): |
| 70 | + data.append( |
| 71 | + { |
| 72 | + "chromosome": chrom, |
| 73 | + "position": positions[i], |
| 74 | + "cumulative_pos": cumulative_positions[i], |
| 75 | + "p_value": pvalues[i], |
| 76 | + "neg_log_p": -np.log10(pvalues[i]), |
| 77 | + } |
| 78 | + ) |
| 79 | + |
| 80 | + cumulative_pos += length * 1e6 |
| 81 | + |
| 82 | +df = pd.DataFrame(data) |
| 83 | + |
| 84 | +# Alternating colors for chromosomes (Python Blue and a gray variant) |
| 85 | +colors = {"odd": "#306998", "even": "#7A9FBF"} |
| 86 | + |
| 87 | +# Create figure |
| 88 | +fig = go.Figure() |
| 89 | + |
| 90 | +# Add scatter traces for each chromosome |
| 91 | +for i, chrom in enumerate(chr_lengths.keys()): |
| 92 | + chr_data = df[df["chromosome"] == chrom] |
| 93 | + color = colors["odd"] if int(chrom) % 2 == 1 else colors["even"] |
| 94 | + |
| 95 | + fig.add_trace( |
| 96 | + go.Scatter( |
| 97 | + x=chr_data["cumulative_pos"], |
| 98 | + y=chr_data["neg_log_p"], |
| 99 | + mode="markers", |
| 100 | + marker=dict(size=5, color=color, opacity=0.7), |
| 101 | + name=f"Chr {chrom}", |
| 102 | + showlegend=False, |
| 103 | + hovertemplate=( |
| 104 | + f"Chr {chrom}<br>Position: %{{customdata[0]:,.0f}} bp<br>-log₁₀(p): %{{y:.2f}}<extra></extra>" |
| 105 | + ), |
| 106 | + customdata=chr_data[["position"]].values, |
| 107 | + ) |
| 108 | + ) |
| 109 | + |
| 110 | +# Genome-wide significance threshold (-log10(5e-8) ≈ 7.3) |
| 111 | +significance_threshold = -np.log10(5e-8) |
| 112 | +fig.add_shape( |
| 113 | + type="line", |
| 114 | + x0=0, |
| 115 | + x1=1, |
| 116 | + xref="paper", |
| 117 | + y0=significance_threshold, |
| 118 | + y1=significance_threshold, |
| 119 | + line=dict(color="#E53935", width=2, dash="dash"), |
| 120 | +) |
| 121 | +fig.add_annotation( |
| 122 | + text="Genome-wide significance (p = 5×10⁻⁸)", |
| 123 | + font=dict(size=16, color="#E53935"), |
| 124 | + xref="paper", |
| 125 | + x=0.99, |
| 126 | + xanchor="right", |
| 127 | + yref="y", |
| 128 | + y=significance_threshold, |
| 129 | + showarrow=False, |
| 130 | + yshift=15, |
| 131 | +) |
| 132 | + |
| 133 | +# Suggestive threshold (-log10(1e-5) = 5) |
| 134 | +suggestive_threshold = 5 |
| 135 | +fig.add_shape( |
| 136 | + type="line", |
| 137 | + x0=0, |
| 138 | + x1=1, |
| 139 | + xref="paper", |
| 140 | + y0=suggestive_threshold, |
| 141 | + y1=suggestive_threshold, |
| 142 | + line=dict(color="#FFD43B", width=2, dash="dot"), |
| 143 | +) |
| 144 | +fig.add_annotation( |
| 145 | + text="Suggestive threshold (p = 10⁻⁵)", |
| 146 | + font=dict(size=16, color="#B8860B"), |
| 147 | + xref="paper", |
| 148 | + x=0.99, |
| 149 | + xanchor="right", |
| 150 | + yref="y", |
| 151 | + y=suggestive_threshold, |
| 152 | + showarrow=False, |
| 153 | + yshift=15, |
| 154 | +) |
| 155 | + |
| 156 | +# Highlight significant SNPs |
| 157 | +significant_snps = df[df["neg_log_p"] > significance_threshold] |
| 158 | +if len(significant_snps) > 0: |
| 159 | + fig.add_trace( |
| 160 | + go.Scatter( |
| 161 | + x=significant_snps["cumulative_pos"], |
| 162 | + y=significant_snps["neg_log_p"], |
| 163 | + mode="markers", |
| 164 | + marker=dict(size=10, color="#E53935", symbol="diamond", line=dict(color="white", width=1)), |
| 165 | + name="Significant SNPs", |
| 166 | + showlegend=True, |
| 167 | + hovertemplate=( |
| 168 | + "Significant SNP<br>" |
| 169 | + "Chr %{customdata[0]}<br>" |
| 170 | + "Position: %{customdata[1]:,.0f} bp<br>" |
| 171 | + "-log₁₀(p): %{y:.2f}<extra></extra>" |
| 172 | + ), |
| 173 | + customdata=significant_snps[["chromosome", "position"]].values, |
| 174 | + ) |
| 175 | + ) |
| 176 | + |
| 177 | +# Chromosome tick positions and labels |
| 178 | +chr_positions = [chr_centers[chrom] for chrom in chr_lengths.keys()] |
| 179 | +chr_labels = list(chr_lengths.keys()) |
| 180 | + |
| 181 | +# Layout |
| 182 | +fig.update_layout( |
| 183 | + title=dict(text="manhattan-gwas · plotly · pyplots.ai", font=dict(size=32), x=0.5, xanchor="center"), |
| 184 | + xaxis=dict( |
| 185 | + title=dict(text="Chromosome", font=dict(size=24)), |
| 186 | + tickfont=dict(size=16), |
| 187 | + tickmode="array", |
| 188 | + tickvals=chr_positions, |
| 189 | + ticktext=chr_labels, |
| 190 | + showgrid=False, |
| 191 | + zeroline=False, |
| 192 | + ), |
| 193 | + yaxis=dict( |
| 194 | + title=dict(text="-log₁₀(p-value)", font=dict(size=24)), |
| 195 | + tickfont=dict(size=18), |
| 196 | + gridcolor="rgba(0,0,0,0.1)", |
| 197 | + gridwidth=1, |
| 198 | + zeroline=False, |
| 199 | + ), |
| 200 | + template="plotly_white", |
| 201 | + legend=dict(yanchor="top", y=0.99, xanchor="left", x=0.01, font=dict(size=16), bgcolor="rgba(255,255,255,0.8)"), |
| 202 | + margin=dict(l=80, r=50, t=80, b=80), |
| 203 | + plot_bgcolor="white", |
| 204 | +) |
| 205 | + |
| 206 | +# Save as PNG (4800 x 2700 px) |
| 207 | +fig.write_image("plot.png", width=1600, height=900, scale=3) |
| 208 | + |
| 209 | +# Save as interactive HTML |
| 210 | +fig.write_html("plot.html", include_plotlyjs="cdn") |
0 commit comments