151 lines
15 KiB
Python
151 lines
15 KiB
Python
"""Post-run descriptive audit only; no selection, simulation or model mutation."""
|
||
from pathlib import Path
|
||
import json,hashlib,importlib.util
|
||
import numpy as np
|
||
import pandas as pd
|
||
|
||
ROOT=Path('/data/digital-drosophila')
|
||
exp=ROOT/'experiments/phase4z-pn-feedback-audit'
|
||
assert not (exp/'checksums.sha256').exists(),'Sealed output is immutable'
|
||
(exp/'audit/descriptive-analysis-provenance.json').write_text(json.dumps({
|
||
'source':str(Path(__file__).resolve()),'sha256':hashlib.sha256(Path(__file__).read_bytes()).hexdigest(),
|
||
'post_run_descriptive_analysis':True,'not_part_of_simulation_preregistration':True,
|
||
'new_simulations_or_selection':False},indent=2)+'\n')
|
||
meta=pd.read_parquet(exp/'audit/incoming-source-metadata.parquet').set_index('neuron_idx')
|
||
watch=pd.read_parquet(exp/'audit/watched-metadata.parquet').set_index('neuron_idx')
|
||
d=pd.read_parquet(exp/'drive/edge-accounting.parquet');d=d[d.period=='late'].copy()
|
||
d['condition']=d.trial.str.split('-',n=2).str[-1]
|
||
rep=pd.read_parquet(exp/'replay/results.parquet');rep['condition']=rep.trial.str.split('-',n=2).str[-1]
|
||
path=pd.read_parquet(exp/'ancestry/event-supported-edges.parquet')
|
||
q=pd.read_parquet(exp/'input-audit/paired-summary.parquet')
|
||
current=json.loads((exp/'protocol/input-hashes.json').read_text())
|
||
prior=json.loads((ROOT/'experiments/phase4y-dpm-hypothesis/protocol/input-hashes.json').read_text())
|
||
old={x['odor']:x for x in prior if x['seed']==131101}
|
||
input_checks=[{'odor_index':x['odor_index'],
|
||
'forcing_hash_exact':x['forcing']==old[x['odor_index']]['forcing_hash'],
|
||
'targets_rates_hash_exact':x['targets_rates']==old[x['odor_index']]['targets_rates_hash']} for x in current]
|
||
assert all(x['forcing_hash_exact'] and x['targets_rates_hash_exact'] for x in input_checks)
|
||
(exp/'regression/input-hashes-vs-Y.json').write_text(json.dumps({'pass':True,'checks':input_checks},indent=2)+'\n')
|
||
rows=[]
|
||
for (trial,node),g in d.groupby(['trial','target_neuron_idx']):
|
||
den=g.positive_arrival_sum.sum()
|
||
for role,z in g.groupby('source_role'):
|
||
rows.append({'trial':trial,'target_neuron_idx':node,'source_role':role,
|
||
'positive_fraction':float(z.positive_arrival_sum.sum()/den) if den else None,
|
||
'positive_arrival_sum':float(z.positive_arrival_sum.sum()),'negative_arrival_sum':float(z.negative_arrival_sum.sum())})
|
||
pd.DataFrame(rows).to_parquet(exp/'drive/late-population-accounting.parquet',index=False)
|
||
top=[];timing=[];transitions=[]
|
||
for trial,g in d[d.target_neuron_idx==82312].groupby('trial'):
|
||
a=g[g.source_role=='ALPN'].sort_values(['positive_arrival_sum','edge_id'],ascending=[False,True]).iloc[0]
|
||
top.append({'trial':trial,'source_neuron_idx':int(a.source_neuron_idx),'body_id':int(meta.loc[a.source_neuron_idx,'body_id']),
|
||
'type':str(meta.loc[a.source_neuron_idx,'annotation_type']),'edge_id':int(a.edge_id),
|
||
'positive_fraction':float(a.positive_arrival_sum/g.positive_arrival_sum.sum()),'edge_event_count':int(a.edge_event_count)})
|
||
pd.DataFrame(top).to_parquet(exp/'drive/dominant-late-PN.parquet',index=False)
|
||
motifs=[]
|
||
for trial,g in path[path.endpoint=='first_late'].groupby('trial'):
|
||
g=g.sort_values('backward_hop')
|
||
pairs=g.groupby(['source_neuron_idx','target_neuron_idx']).size()
|
||
for (src,dst),n in pairs.items():
|
||
if (dst,src) in pairs.index:
|
||
motifs.append({'trial':trial,'source_neuron_idx':int(src),'target_neuron_idx':int(dst),
|
||
'forward_occurrences_on_ancestry':int(n),'reverse_occurrences_on_ancestry':int(pairs.loc[(dst,src)]),
|
||
'claim':'actual time-supported reciprocal motif on ranked path, not unique necessity'})
|
||
pd.DataFrame(motifs).to_parquet(exp/'ancestry/reciprocal-event-motifs.parquet',index=False)
|
||
loop=path[(path.endpoint=='first_late')&path.source_neuron_idx.isin([1381,1511])&path.target_neuron_idx.isin([1381,1511])]
|
||
for trial,g in loop.groupby('trial'):
|
||
timing.append({'trial':trial,'recurrent_edge_events_on_path':len(g),
|
||
'first_emission_ms':float(g.source_step.min()/10),'last_arrival_ms':float(g.arrival_step.max()/10),
|
||
'all_delays_exact':bool(((g.arrival_step-g.source_step)==18).all()),
|
||
'all_arrivals_precede_target':bool((g.arrival_step<=g.target_step).all())})
|
||
pd.DataFrame(timing).to_parquet(exp/'ancestry/loop-timing.parquet',index=False)
|
||
for r in q.itertuples():
|
||
z=pd.read_parquet(exp/'input-audit'/f'{r.odor_index}-{r.seed}-paired-requests.parquet')
|
||
diff=z[z.status_changed]
|
||
for (before,after),g in diff.groupby(['status_reference','status_off']):
|
||
transitions.append({'odor_index':r.odor_index,'reference_status':before,'off_status':after,'count':len(g)})
|
||
pd.DataFrame(transitions).to_parquet(exp/'input-audit/status-transitions.parquet',index=False)
|
||
isirows=[]
|
||
for f in sorted((exp/'states').glob('*-watched.npz')):
|
||
trial=f.name.removesuffix('-watched.npz');z=np.load(f);st=z['state'];ids=z['watch']
|
||
for j,node in enumerate(ids):
|
||
steps=np.flatnonzero(st[5000:,6,j])+5000;isi=np.diff(steps)
|
||
isirows.append({'trial':trial,'neuron_idx':int(node),'late_natural_spikes':len(steps),
|
||
'late_rate_hz':2*len(steps),'ISI_min_ms':float(isi.min()/10) if len(isi) else None,
|
||
'ISI_mean_ms':float(isi.mean()/10) if len(isi) else None,'ISI_max_ms':float(isi.max()/10) if len(isi) else None,
|
||
'ISI_equal_refractory_fraction':float(np.mean(isi==22)) if len(isi) else None,
|
||
'late_post_g_mean':float(st[5000:,5,j].mean(dtype=np.float64)),
|
||
'late_post_g_min':float(st[5000:,5,j].min()),'late_post_g_max':float(st[5000:,5,j].max())})
|
||
isidf=pd.DataFrame(isirows);isidf.to_parquet(exp/'drive/late-ISI-and-g.parquet',index=False)
|
||
# Read-only coverage audit of the frozen graded APL replacement, not an intervention.
|
||
canonical=ROOT/'datasets/male-cns-v1.0/processed/neuron-only'
|
||
fullmeta=pd.read_parquet(canonical/'neuron-index/unified-metadata.parquet').sort_values('neuron_idx').reset_index(drop=True)
|
||
assert np.array_equal(fullmeta.neuron_idx.to_numpy(),np.arange(166700))
|
||
ptr=np.load(canonical/'gpu/csr-indptr.int64.npy',mmap_mode='r')
|
||
dst=np.load(canonical/'gpu/csr-indices.int32.npy',mmap_mode='r')
|
||
raw=np.load(canonical/'gpu/csr-weights.int64.npy',mmap_mode='r')
|
||
aplids=np.flatnonzero(fullmeta.annotation_type.eq('APL'));assert len(aplids)==2
|
||
aplrows=[]
|
||
for src in aplids:
|
||
for eid in range(int(ptr[src]),int(ptr[src+1])):
|
||
node=int(dst[eid]);m=fullmeta.iloc[node]
|
||
aplrows.append({'edge_id':eid,'source_neuron_idx':int(src),'source_body_id':int(fullmeta.iloc[src].body_id),
|
||
'target_neuron_idx':node,'target_body_id':int(m.body_id),'target_class':m.annotation_class,
|
||
'target_type':m.annotation_type,'raw_synapse_count':int(raw[eid]),
|
||
'frozen_APL_fast_edge_zeroed':True,'graded_replacement_target':bool(m.annotation_class=='Kenyon_Cell')})
|
||
apldf=pd.DataFrame(aplrows);apldf.to_parquet(exp/'audit/APL-target-coverage.parquet',index=False)
|
||
omitted=apldf[~apldf.graded_replacement_target]
|
||
focus=omitted[omitted.target_neuron_idx.isin(watch[watch.annotation_class=='ALPN'].index)]
|
||
spec=importlib.util.find_spec('phase4p_revalidate');assert spec and spec.origin
|
||
srcfile=Path(spec.origin);assert srcfile.exists()
|
||
aplsummary={'scope':'static coverage of frozen DeterministicAPLSim.install and graded output matrix; not proof of functional effect',
|
||
'APL_neurons':len(aplids),'all_outgoing_edges':len(apldf),'all_outgoing_synapse_count':int(apldf.raw_synapse_count.sum()),
|
||
'non_KC_edges_without_graded_replacement':len(omitted),'non_KC_synapse_count_without_graded_replacement':int(omitted.raw_synapse_count.sum()),
|
||
'non_KC_target_neurons':int(omitted.target_neuron_idx.nunique()),'watched_PN_edges':len(focus),
|
||
'watched_PN_synapse_count':int(focus.raw_synapse_count.sum()),'source':str(srcfile),
|
||
'source_sha256':hashlib.sha256(srcfile.read_bytes()).hexdigest(),'baseline_not_modified':True,
|
||
'cannot_replace_terminal_inhibition_with_whole_PN_inhibition_without_new_hypothesis':True}
|
||
(exp/'audit/APL-target-coverage.json').write_text(json.dumps(aplsummary,indent=2)+'\n')
|
||
lines=['# Phase 4Z — mechanism interpretation','',
|
||
'Post-run descriptive analysis, not a new candidate gate. No tuning, learning, new stimulus or sign hypothesis.',
|
||
'The earlier generic label PN feedback is refined below using actual ALLN cell identities. Do not assume the residual loop is KC→MBON→PN feedback merely because that pathway exists in the literature.','',
|
||
'## Watched identities','', '| neuron_idx | body_id | class | type | metadata transmitter | raw ground_truth |','|---:|---:|---|---|---|---|']
|
||
for idx,r in watch.iterrows():lines.append(f'| {idx} | {r.body_id} | {r.annotation_class} | {r.annotation_type} | {r.neurotransmitter} | {r.nt_ground_truth} |')
|
||
lines+=['','## Representative first-late ancestry','',
|
||
'Odor 0 / seed 131101 / frozen reference. Ordered source spike → edge arrival → target spike (ms). These are ranked event-supported contributions; no unique predecessor assertion.',
|
||
'| source idx/type | target idx/type | source time | arrival time | target time | contacts | effective weight |','|---|---|---:|---:|---:|---:|---:|']
|
||
example=path[(path.trial=='0-131101-frozen_reference')&(path.endpoint=='first_late')].sort_values('backward_hop').head(5).iloc[::-1]
|
||
for r in example.itertuples():lines.append(f'| {r.source_neuron_idx}/{r.source_type} | {r.target_neuron_idx}/{r.target_type} | {r.source_step/10:.1f} | {r.arrival_step/10:.1f} | {r.target_step/10:.1f} | {r.raw_synapse_count} | {r.effective_weight:.6f} |')
|
||
lines+=['','## Dominant PN arrival contribution','', '| trial | PN idx/type | positive arrival fraction |','|---|---|---:|']
|
||
for r in top:lines.append(f'| {r["trial"]} | {r["source_neuron_idx"]}/{r["type"]} | {100*r["positive_fraction"]:.3f}% |')
|
||
lines+=['','Fractions are sums of positive synaptic arrivals, not percentages of causal responsibility. Negative input, recorded APL and membrane dynamics are retained in exact CPU replay.',
|
||
f'The preidentified 1381↔1511 event motif occurs in {len(timing)}/12 first-late ancestry paths. Full temporal edges and alternative inputs are retained.',
|
||
'Conditional omission results are in replay/results.parquet and REPORT.md. They retain other source events; loss of spikes there is not evidence that a whole-brain ablation would produce the same result.',
|
||
'Critically, omitting ALLN arrivals to a watched ALLN does not silence it in conditional replay: recorded PN and other inputs can still sustain high firing. Thus the 1381↔1511 motif cannot be called the sole cause or an isolated autonomous oscillator.',
|
||
'Removing the existing VA6_adPN→KC82312 edge in conditional replay eliminates target late firing; removing ALLN input to VA6_adPN also eliminates its late firing. These are two conditional links, not a demonstrated whole-brain necessity chain.','',
|
||
'## Late firing of preidentified ALLNs','', '| neuron_idx | late rate range Hz | min/max ISI ms | fraction ISI = refractory, range |','|---:|---:|---:|---:|']
|
||
for node,g in isidf[isidf.neuron_idx.isin([1197,1381,1511])].groupby('neuron_idx'):
|
||
lines.append(f'| {node} | {g.late_rate_hz.min()}–{g.late_rate_hz.max()} | {g.ISI_min_ms.min():.2f}/{g.ISI_max_ms.max():.2f} | {g.ISI_equal_refractory_fraction.min():.5f}–{g.ISI_equal_refractory_fraction.max():.5f} |')
|
||
lines+=['','This is a measured property of the mathematical model, not a validated biological firing rate. Window-boundary counting can put spikes/window slightly above 1/refractory. No new biological runaway threshold is inferred.','',
|
||
'## Frozen APL replacement target coverage','',
|
||
f'The frozen implementation zeroes all {len(apldf)} APL fast outgoing edges and reintroduces graded inhibition only on KC targets. {len(omitted)} non-KC edges ({int(omitted.raw_synapse_count.sum())} anatomical contacts) have no direct graded replacement in this implementation.',
|
||
f'This includes {len(focus)} edges / {int(focus.raw_synapse_count.sum())} anatomical contacts onto the seven watched PNs. Raw edges and counts remain intact. These edges were not removed by this audit.',
|
||
'This is an implementation coverage gap relative to the anatomical graph, not evidence that restoring a whole-neuron negative current would be biologically correct. APL inhibition of PN boutons is spatial/terminal-specific and a point-neuron approximation may not capture it.',
|
||
'Primary study: Prisco et al., eLife 10:e74172 (2021), https://pubmed.ncbi.nlm.nih.gov/34964714/ . It reports reciprocal APL inhibition targeting PN boutons and KC dendrites; it does not provide these MaleCNS edge-specific effective parameters.',
|
||
'No frozen APL parameter, edge, sign or output target was modified. Functional necessity of this coverage gap remains untested.','',
|
||
'## What input matching means','',
|
||
f'All six paired request schedules match. Status differs for {int(q.status_changed_requests.sum())} requests in total; actual output at requested times differs for {int(q.delivered_output_differences.sum())} requests.',
|
||
'An accepted→preexisting transition changes the accepted-new count without removing the output spike. A blocked transition can change actual delivery. Thus accepted-new totals alone are insufficient to describe input matching.',
|
||
'The first differing request is state dependent; this is an endogenous response to an existing intervention, not evidence of different external RNG. No clamp, injected substitute events or refractory override was introduced.',
|
||
'The earlier closed-loop intervention remains interpretable as a total system effect under matched requests. It is not an estimate of a direct DPM effect holding actual PN output trains fixed.','',
|
||
'## Literature boundary','',
|
||
'Excitatory local antennal-lobe circuits and reciprocal local-neuron/PN excitation have experimental support. They must not automatically be relabeled inhibitory because the model has persistent activity.',
|
||
'Shang et al., Cell 2007: https://pubmed.ncbi.nlm.nih.gov/17289577/',
|
||
'Huang et al., Neuron 2010, functional connectivity of excitatory local interneurons: https://doi.org/10.1016/j.neuron.2010.08.025',
|
||
'These experiments do not validate uniform current weights, this MaleCNS individual, or its simulated persistent firing. No receptor, adaptation or electrical-synapse parameters are inferred or inserted.',
|
||
'Next priority is independent assessment of local AL-loop dynamics and input-interface causal controls, not changing dopamine learning rate. No official model replacement is justified.']
|
||
(exp/'audit/MECHANISM.md').write_text('\n'.join(lines)+'\n',encoding='utf-8')
|
||
print(json.dumps({'dominant_PNs':sorted({r['source_neuron_idx'] for r in top}),
|
||
'dominant_PN_positive_fraction_range':[min(r['positive_fraction'] for r in top),max(r['positive_fraction'] for r in top)],
|
||
'trials_with_preidentified_reciprocal_motif':len(timing),
|
||
'status_changed_requests':int(q.status_changed_requests.sum()),'delivered_output_differences':int(q.delivered_output_differences.sum())}),flush=True)
|