import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.gridspec as gridspec
import seaborn as sns
from pathlib import Path

pd.set_option('display.max_colwidth', 90)
pd.set_option('display.width', 200)
sns.set_context('notebook')
sns.set_style('whitegrid')

PROJ = Path('/home/aparkin/BERIL-research-observatory/projects/caulobacter_fur_lipida_loss')
DATA = PROJ / 'data'
FIG = PROJ / 'figures'

# Load saved outputs
nb01 = pd.read_csv(DATA / 'NB01_fur_only_signature.csv')
nb02_A = pd.read_csv(DATA / 'NB02_pathA_concordant_strong_scoring.csv')
nb02_B = pd.read_csv(DATA / 'NB02_pathB_buffered_scoring.csv')
nb03_cohorts = pd.read_csv(DATA / 'NB03_chvi_phase_cohorts.csv')
nb03_themes = pd.read_csv(DATA / 'NB03_phase_cohort_themes.csv')
nb04_sph = pd.read_csv(DATA / 'NB04_sphingolipid_transcript.csv')
nb04_lpt = pd.read_csv(DATA / 'NB04_lpt_transcript.csv')
nb04_prot = pd.read_csv(DATA / 'NB04_protein_panel.csv')
nb05_hits = pd.read_csv(DATA / 'NB05_pg_significant_hits.csv')
nb06_bool = pd.read_csv(DATA / 'NB06_comparative_presence_bool.csv', index_col=0)
nb06_counts = pd.read_csv(DATA / 'NB06_comparative_presence_counts.csv')

print('Loaded all NB01-NB06 outputs.')

Loaded all NB01-NB06 outputs.

# Use a tight grid layout
fig = plt.figure(figsize=(18, 14))
gs = gridspec.GridSpec(2, 2, figure=fig, hspace=0.35, wspace=0.32)

# ========== PANEL A: Δfur signature concordance ==========
axA = fig.add_subplot(gs[0, 0])

both = nb01.dropna(subset=['log2fc_ours','log2fc_fur']).copy()
# Categorize
def cat(row):
    if pd.isna(row['log2fc_ours']): return 'missing'
    same = np.sign(row['log2fc_fur']) == np.sign(row['log2fc_ours'])
    ao = abs(row['log2fc_ours'])
    if not same and ao > 0.5: return 'discordant'
    if same and ao > 1: return 'concordant_strong'
    if same and ao > 0.5: return 'concordant_weak'
    return 'buffered'
both['cat'] = both.apply(cat, axis=1)

color_map = {'concordant_strong':'#1f77b4', 'concordant_weak':'#aec7e8',
             'buffered':'#ff7f0e', 'discordant':'#d62728', 'missing':'lightgrey'}
for cat_, sub in both.groupby('cat'):
    axA.scatter(sub['log2fc_fur'], sub['log2fc_ours'],
                color=color_map.get(cat_, 'grey'),
                label=f'{cat_} (n={len(sub)})', s=45, alpha=0.8,
                edgecolor='k', linewidth=0.3)

lim = 11
axA.plot([-lim, lim], [-lim, lim], '--', color='red', alpha=0.4, lw=1)
axA.axhline(0, color='grey', lw=0.5); axA.axvline(0, color='grey', lw=0.5)
axA.set_xlim(-lim, lim); axA.set_ylim(-lim, lim)
axA.set_xlabel('Leaden 2018 — log2FC Δfur / WT', fontsize=11)
axA.set_ylabel('Our data — log2FC 4584 / 4580\n(Δfur ΔsspB ΔrsaA vs ΔrsaA)', fontsize=11)
axA.set_title('(A) Δfur signature: concordance and ΔsspB-buffered set\nSpearman ρ = 0.315, p = 2.08e-03 over 93 Leaden Δfur DEGs',
              fontsize=12, fontweight='bold')
axA.legend(fontsize=9, loc='lower right', framealpha=0.9)

# Annotate key genes
key_anns = {'CCNA_02277':'HutA', 'CCNA_00055':'fur', 'CCNA_01475':'OmpW', 'CCNA_01467':'ccoN'}
for ccna, label in key_anns.items():
    row = both[both['locustag']==ccna]
    if len(row):
        axA.annotate(label, (row['log2fc_fur'].iloc[0], row['log2fc_ours'].iloc[0]),
                    xytext=(5,5), textcoords='offset points', fontsize=9, fontweight='bold')


# ========== PANEL B: ChvI phase cohort theme shift ==========
axB = fig.add_subplot(gs[0, 1])

# Bar chart: themes per phase cohort (early vs late)
themes_subset = nb03_themes[~nb03_themes['theme'].isin(['hypothetical'])].copy()
themes_long = pd.DataFrame({
    'theme': list(themes_subset['theme']) * 2,
    'cohort': ['early'] * len(themes_subset) + ['late'] * len(themes_subset),
    'pct': list(themes_subset['early_pct']) + list(themes_subset['late_pct'])
})

theme_order = themes_subset.sort_values('late_pct', ascending=True)['theme'].tolist()
themes_long['theme'] = pd.Categorical(themes_long['theme'], categories=theme_order, ordered=True)

palette = {'early':'#74a9cf', 'late':'#08519c'}
sns.barplot(data=themes_long, y='theme', x='pct', hue='cohort', palette=palette, ax=axB)
axB.set_xlabel('% of cohort with theme annotation', fontsize=11)
axB.set_ylabel('')
axB.set_title('(B) ChvI regulon phase shift\n(early = 30 genes induced in 4584-vs-4580; late = 49 genes induced only in 4599-vs-4584)',
              fontsize=12, fontweight='bold')
axB.legend(title='Phase', loc='lower right')
axB.grid(True, axis='x', alpha=0.3)


# ========== PANEL C: Sphingolipid + Lpt + CtpA ==========
axC = fig.add_subplot(gs[1, 0])

# Build summary points: each gene with logFC in 4584-vs-4580 and 4599-vs-4584
def build_summary(df, set_label, color, marker, role_col='role'):
    df = df.copy()
    df['set'] = set_label
    df['color'] = color
    df['marker'] = marker
    return df[['locustag','role','set','color','marker','logFC_4584_vs_4580','logFC_4599_vs_4584']]

# Pick representatives for visual clarity
sph_show = ['CCNA_01220','CCNA_01218','CCNA_01212','CCNA_01213','CCNA_01214','CCNA_01226']
lpt_show = ['CCNA_00307','CCNA_01760','CCNA_03716']
ctpa_show = ['CCNA_03113']

sph_sub = nb04_sph[nb04_sph['locustag'].isin(sph_show)].assign(set_label='sphingolipid', color='#1b9e77')
lpt_sub = nb04_lpt[nb04_lpt['locustag'].isin(lpt_show)].assign(set_label='Lpt apparatus', color='#7570b3')
ctpa_sub = pd.read_csv(DATA / 'NB04_ctpA_transcript.csv').assign(set_label='CtpA', color='#d95f02')

for sub in [sph_sub, lpt_sub, ctpa_sub]:
    axC.scatter(sub['logFC_4584_vs_4580'], sub['logFC_4599_vs_4584'],
                color=sub['color'].iloc[0], s=120, alpha=0.85,
                edgecolor='k', linewidth=0.5,
                label=f'{sub["set_label"].iloc[0]} (n={len(sub)})')
    for _, row in sub.iterrows():
        # Short annotation
        label = row['role'].split('—')[0].strip()[:12]
        axC.annotate(label, (row['logFC_4584_vs_4580'], row['logFC_4599_vs_4584']),
                    xytext=(5,3), textcoords='offset points', fontsize=8)

axC.axhline(0, color='grey', lw=0.5); axC.axvline(0, color='grey', lw=0.5)
axC.set_xlabel('log2FC 4584-vs-4580  (Δfur ΔsspB effect)', fontsize=11)
axC.set_ylabel('log2FC 4599-vs-4584  (added Δlpxc effect)', fontsize=11)
axC.set_title('(C) Sphingolipid + Lpt + CtpA panel\nSphingolipid pathway constitutive; Lpt apparatus maintained or up; lptC2 protein +2.1x',
              fontsize=12, fontweight='bold')
axC.legend(loc='upper left', fontsize=10, framealpha=0.9)
axC.set_xlim(-1.2, 1.2); axC.set_ylim(-1.2, 1.2)


# ========== PANEL D: Comparative species presence/absence ==========
axD = fig.add_subplot(gs[1, 1])

# Select key families for the panel — focus on storytelling rather than completeness
key_families = [
    'spt — serine palmitoyltransferase',
    'bcerS — bacterial ceramide synthase',
    'sphingosine kinase (sphk)',
    'ChvG / ExoS — sensor kinase',
    'ChvI / ExoR — response regulator',
    'Capsule biosynthesis (siaD/cps)',
    'Polysialyltransferase',
    'Late acyltransferase (lpxX/lpxL)',
    'PBP1A / MrcA — class A penicillin-binding',
    'Ld-transpeptidase (LdtJ/K-class)',
]
panel_d = nb06_bool.loc[[f for f in key_families if f in nb06_bool.index]]
panel_d.columns = [c.replace('_present','') for c in panel_d.columns]

sns.heatmap(panel_d, cmap=['#fed976','#08519c'], cbar=False, annot=panel_d.values, fmt='d',
            linewidths=0.5, linecolor='white',
            yticklabels=panel_d.index, xticklabels=panel_d.columns, ax=axD)
axD.set_title('(D) Comparative presence/absence (PaperBLAST)\nSphingolipid biosynthesis + ChvI are uniquely Caulobacter',
              fontsize=12, fontweight='bold')
axD.set_xlabel(''); axD.set_ylabel('')
for label in axD.get_yticklabels():
    label.set_fontsize(9)

fig.suptitle('Caulobacter Δfur Δsspb permits Δlpxc — multi-level mechanism (Phase A finding for Project caulobacter_fur_lipida_loss)',
             fontsize=14, fontweight='bold', y=0.995)
plt.savefig(FIG / 'NB07_synthesis_master.png', dpi=160, bbox_inches='tight')
plt.show()

scorecard = pd.DataFrame([
    {'Hypothesis':'H1 — ChvI cooperator + consequence',
     'Sub-claim':'Phase structure (≥10 per cohort)',
     'Threshold':'early ≥10 AND late ≥10',
     'Observed':'early=30, late=49',
     'Verdict':'PASS'},
    {'Hypothesis':'H1 — ChvI cooperator + consequence',
     'Sub-claim':'SigU drives late cohort',
     'Threshold':'late ≥50% envelope/transport/regulator',
     'Observed':'24.5% (literature gap)',
     'Verdict':'PARTIAL'},
    {'Hypothesis':'H2 — Critical Fur regulon subset',
     'Sub-claim':'Path A (concordant_strong)',
     'Threshold':'≥10% phenotype-bearing',
     'Observed':'17/32 = 53%',
     'Verdict':'PASS strongly'},
    {'Hypothesis':'H2 — Critical Fur regulon subset',
     'Sub-claim':'Path B (SspB-buffered cbb3/fix)',
     'Threshold':'≥10% phenotype-bearing',
     'Observed':'9/26 = 35%',
     'Verdict':'PASS'},
    {'Hypothesis':'H3 — Sphingolipid + Lpt repurposing',
     'Sub-claim':'CtpA upregulation',
     'Threshold':'FDR<0.05 in 4599-vs-4584 OR ≥2x protein',
     'Observed':'p=0.11 transcript; not in OM proteome; cumulative 4599-vs-4580 FDR=0.035',
     'Verdict':'BORDERLINE'},
    {'Hypothesis':'H3 — Sphingolipid + Lpt repurposing',
     'Sub-claim':'Sphingolipid pathway constitutive',
     'Threshold':'No biosynthesis gene logFC>1, FDR<0.05',
     'Observed':'0 up, several mildly down (spt -0.64, sphk -0.40)',
     'Verdict':'PASS strongly'},
    {'Hypothesis':'H3 — Sphingolipid + Lpt repurposing',
     'Sub-claim':'Canonical Lpt apparatus maintained',
     'Threshold':'No component logFC<-0.5, FDR<0.05',
     'Observed':'0 down; MsbA-like +0.89 FDR=0.01, LptC-related +0.56 FDR=0.005',
     'Verdict':'PASS strongly'},
    {'Hypothesis':'H3 — Sphingolipid + Lpt repurposing',
     'Sub-claim':'(bonus) lptC2 protein-level induction',
     'Threshold':'(post-hoc finding)',
     'Observed':'transcript -0.60 (FDR 0.034), protein log2 +1.08 (>2x)',
     'Verdict':'NOVEL FINDING'},
    {'Hypothesis':'H4 — PG remodeling participates',
     'Sub-claim':'PG enzymes engaged in 4599-vs-4584',
     'Threshold':'≥3 enzymes FDR<0.05 OR |log2|>1 protein',
     'Observed':'31 total (26 transcript, 7 protein)',
     'Verdict':'PASS strongly'},
    {'Hypothesis':'H4 — PG remodeling participates',
     'Sub-claim':'(novel) Pal-Tol envelope-integrity engaged',
     'Threshold':'(post-hoc finding)',
     'Observed':'Pal CCNA_00784 +2.08 transcript, +2.84 protein',
     'Verdict':'NOVEL FINDING'},
    {'Hypothesis':'Comparative species arm',
     'Sub-claim':'Sphingolipid pathway Caulobacter-unique',
     'Threshold':'(presence/absence)',
     'Observed':'spt, bcerS, sphk absent in A.b., N.m., M.c.',
     'Verdict':'CONFIRMED'},
    {'Hypothesis':'Comparative species arm',
     'Sub-claim':'ChvG-ChvI alphaproteobacterial only',
     'Threshold':'(presence/absence)',
     'Observed':'absent in A.b., N.m., M.c.',
     'Verdict':'CONFIRMED'},
])

print('=== HYPOTHESIS SCORECARD ===')
display(scorecard)
scorecard.to_csv(DATA / 'NB07_scorecard.csv', index=False)

=== HYPOTHESIS SCORECARD ===

print('=== NB07 SYNTHESIS COMPLETE ===\n')
print('Master figure: figures/NB07_synthesis_master.png')
print('Scorecard:     data/NB07_scorecard.csv')
print()
print('All 7 notebooks executed; H1-H4 + comparative arm covered.')
print('Ready for /synthesize to write REPORT.md.')

=== NB07 SYNTHESIS COMPLETE ===

Master figure: figures/NB07_synthesis_master.png
Scorecard:     data/NB07_scorecard.csv

All 7 notebooks executed; H1-H4 + comparative arm covered.
Ready for /synthesize to write REPORT.md.

07 Synthesis

07 — Synthesis: master figure summarizing the rescue mechanism¶

Purpose¶

Panels¶

Build the 4-panel master figure¶

Hypothesis scorecard¶

Mechanism summary¶

	Hypothesis	Sub-claim	Threshold	Observed	Verdict
0	H1 — ChvI cooperator + consequence	Phase structure (≥10 per cohort)	early ≥10 AND late ≥10	early=30, late=49	PASS
1	H1 — ChvI cooperator + consequence	SigU drives late cohort	late ≥50% envelope/transport/regulator	24.5% (literature gap)	PARTIAL
2	H2 — Critical Fur regulon subset	Path A (concordant_strong)	≥10% phenotype-bearing	17/32 = 53%	PASS strongly
3	H2 — Critical Fur regulon subset	Path B (SspB-buffered cbb3/fix)	≥10% phenotype-bearing	9/26 = 35%	PASS
4	H3 — Sphingolipid + Lpt repurposing	CtpA upregulation	FDR<0.05 in 4599-vs-4584 OR ≥2x protein	p=0.11 transcript; not in OM proteome; cumulative 4599-vs-4580 FDR=0.035	BORDERLINE
5	H3 — Sphingolipid + Lpt repurposing	Sphingolipid pathway constitutive	No biosynthesis gene logFC>1, FDR<0.05	0 up, several mildly down (spt -0.64, sphk -0.40)	PASS strongly
6	H3 — Sphingolipid + Lpt repurposing	Canonical Lpt apparatus maintained	No component logFC<-0.5, FDR<0.05	0 down; MsbA-like +0.89 FDR=0.01, LptC-related +0.56 FDR=0.005	PASS strongly
7	H3 — Sphingolipid + Lpt repurposing	(bonus) lptC2 protein-level induction	(post-hoc finding)	transcript -0.60 (FDR 0.034), protein log2 +1.08 (>2x)	NOVEL FINDING
8	H4 — PG remodeling participates	PG enzymes engaged in 4599-vs-4584	≥3 enzymes FDR<0.05 OR \|log2\|>1 protein	31 total (26 transcript, 7 protein)	PASS strongly
9	H4 — PG remodeling participates	(novel) Pal-Tol envelope-integrity engaged	(post-hoc finding)	Pal CCNA_00784 +2.08 transcript, +2.84 protein	NOVEL FINDING
10	Comparative species arm	Sphingolipid pathway Caulobacter-unique	(presence/absence)	spt, bcerS, sphk absent in A.b., N.m., M.c.	CONFIRMED
11	Comparative species arm	ChvG-ChvI alphaproteobacterial only	(presence/absence)	absent in A.b., N.m., M.c.	CONFIRMED