1import heapq
  2import json
  3import time
  4import numpy as np
  5from sklearn.datasets import load_digits
  6from sklearn.model_selection import train_test_split
  7from sklearn.linear_model import LogisticRegression
  8from sklearn.metrics import accuracy_score
  9
 10SEED = 7
 11
 12class EventGraph:
 13    def __init__(self, n_source, n_hidden, n_out, source_edges, hidden_edges, vmax=100):
 14        self.S, self.H, self.O = n_source, n_hidden, n_out
 15        self.source_edges = source_edges
 16        self.hidden_edges = hidden_edges
 17        self.vmax = vmax
 18
 19    def run(self, x, return_trace=False):
 20        a = np.zeros(self.H, dtype=float)
 21        visits = np.zeros(self.H, dtype=int)
 22        z = np.zeros(self.O, dtype=float)
 23        q, serial = [], 0
 24        for i, value in enumerate(x):
 25            heapq.heappush(q, (0, serial, i, float(value)))
 26            serial += 1
 27        events = 0
 28        trace = []
 29        while q and events < 100000:
 30            t, _, node, value = heapq.heappop(q)
 31            events += 1
 32            if node < self.S:
 33                for si, h, w, d in self.source_edges:
 34                    if si == node:
 35                        heapq.heappush(q, (t + d, serial, self.S + h, w * value))
 36                        serial += 1
 37                continue
 38            h = node - self.S
 39            if h >= self.H:
 40                z[h - self.H] += value
 41                continue
 42            if visits[h] >= self.vmax:
 43                continue
 44            a[h] += value
 45            visits[h] += 1
 46            y = np.tanh(a[h])
 47            trace.append((t, h, a[h], y, visits[h]))
 48            for hh, out, w, d in self.hidden_edges:
 49                if hh == h:
 50                    heapq.heappush(q, (t + d, serial, self.S + self.H + out, w * y))
 51                    serial += 1
 52        result = dict(z=z, accum=a, visits=visits, events=events, trace=trace)
 53        return result if return_trace else z
 54
 55
 56def prediction_sweeps():
 57    rng = np.random.default_rng(SEED)
 58    rows = []
 59    # Prediction 1: with one shared node and linear response, final accumulator is sum of arrivals.
 60    # For k equal arrivals c, a = k*c exactly; measured slope should be 1.
 61    ks = np.arange(1, 9)
 62    measured = []
 63    for k in ks:
 64        edges = [(i, 0, 1.0, int(i % 2)) for i in range(k)]
 65        g = EventGraph(k, 1, 1, edges, [(0, 0, 1.0, 0)])
 66        measured.append(g.run(np.full(k, 0.1), True)['accum'][0])
 67    slope = float(np.polyfit(ks, measured, 1)[0] / 0.1)
 68    rows.append({'prediction': 'shared accumulator slope versus number of equal arrivals = 1/cancel-free sum',
 69                  'predicted': 1.0, 'observed': slope, 'relative_error': abs(slope-1.0)})
 70
 71    # Prediction 2: tanh response saturates; equal arrivals k*c approach 1, with
 72    # inverse-tanh(0.9)/c predicted for the 0.9 crossing.
 73    c = 0.2
 74    crossing = None
 75    ys = []
 76    for k in range(1, 31):
 77        edges = [(i, 0, 1.0, 0) for i in range(k)]
 78        g = EventGraph(k, 1, 1, edges, [(0, 0, 1.0, 0)])
 79        y = g.run(np.full(k, c), True)['trace'][-1][3]
 80        ys.append(y)
 81        if crossing is None and y >= .9:
 82            crossing = k
 83    predicted_crossing = int(np.ceil(np.arctanh(.9) / c))
 84    rows.append({'prediction': 'shared tanh response reaches 0.9 at k=ceil(atanh(0.9)/c)',
 85                  'predicted': predicted_crossing, 'observed': crossing,
 86                  'relative_error': abs(crossing-predicted_crossing)/predicted_crossing})
 87
 88    # Prediction 3: visit cap Vmax truncates exactly after Vmax hidden visits.
 89    cap_rows = []
 90    k = 8
 91    for cap in range(1, 6):
 92        edges = [(i, 0, 1.0, 0) for i in range(k)]
 93        g = EventGraph(k, 1, 1, edges, [(0, 0, 1.0, 0)], vmax=cap)
 94        r = g.run(np.ones(k), True)
 95        cap_rows.append((cap, int(r['visits'][0]), float(r['accum'][0])))
 96    rows.append({'prediction': 'visit cap produces exactly min(k,Vmax) accepted updates',
 97                 'predicted': [min(k, c) for c in range(1, 6)],
 98                 'observed': [x[1] for x in cap_rows],
 99                 'relative_error': 0.0 if all(x[1] == min(k, x[0]) for x in cap_rows) else 1.0})
100    return rows
101
102
103def make_features(X, shared=True, P=None, Q=None, seed=SEED):
104    rng = np.random.default_rng(seed)
105    n, d = X.shape
106    S, H, O = 16, 24, 10
107    if P is None:
108        P = rng.normal(0, 1/np.sqrt(d), (d, S))
109    if Q is None:
110        Q = rng.normal(0, 1/np.sqrt(H), (H, O))
111    out = []
112    event_counts = []
113    for x in X:
114        src = x @ P
115        if shared:
116            edges = [(i, int(i % 8), 0.55, i % 3) for i in range(S)]
117        else:
118            edges = [(i, i, 0.55, 0) for i in range(S)]
119        hidden_edges = [(h, o, Q[h, o], 0) for h in range(H) for o in range(O)]
120        g = EventGraph(S, H, O, edges, hidden_edges, vmax=20)
121        r = g.run(src, return_trace=True)
122        out.append(r['z'])
123        event_counts.append(r['events'])
124    return np.asarray(out), float(np.mean(event_counts))
125
126
127def classification_test():
128    digits = load_digits()
129    X = digits.data.astype(float) / 16.0
130    y = digits.target
131    Xtr, Xte, ytr, yte = train_test_split(X, y, test_size=.25, random_state=SEED, stratify=y)
132    rng = np.random.default_rng(SEED)
133    P = rng.normal(0, 1/np.sqrt(X.shape[1]), (X.shape[1], 16))
134    Q = rng.normal(0, 1/np.sqrt(24), (24, 10))
135    result = {}
136    for shared in (False, True):
137        t0 = time.perf_counter()
138        Ftr, etr = make_features(Xtr, shared=shared, P=P, Q=Q)
139        Fte, ete = make_features(Xte, shared=shared, P=P, Q=Q)
140        clf = LogisticRegression(max_iter=300, C=1.0, random_state=SEED)
141        clf.fit(Ftr, ytr)
142        pred = clf.predict(Fte)
143        result['shared' if shared else 'no_sharing'] = {
144            'accuracy': float(accuracy_score(yte, pred)),
145            'seconds': time.perf_counter() - t0,
146            'train_mean_events': etr, 'test_mean_events': ete,
147            'feature_dim': int(Ftr.shape[1])
148        }
149    return result
150
151
152def main():
153    report = {'seed': SEED, 'mechanism_checks': prediction_sweeps(), 'classification': classification_test()}
154    with open('results.json', 'w') as f:
155        json.dump(report, f, indent=2)
156    print(json.dumps(report, indent=2))
157
158if __name__ == '__main__':
159    main()