Retina Laboratory / retina-unified-1.0.0

Models &
verification.

Equations, assumptions and independent numerical evidence.

Unified experimental scope

The Retina Laboratory (retina-unified-1.0.0) uses one scene, retinal coordinate chart, sampled cell population and experiment definition for electrical, PRIMA photovoltaic, direct optogenetic and electrical–optogenetic hybrid stimulation. The new biological coupling is a research prototype, not a patient-calibrated prediction. Numerically verified fields and circuits do not validate the neurons or the appearance of human vision.

The retained methods remain available:

Common coordinates and anatomy

All distances are in µm; time is ms. The fovea is (0,0), retinal x points nasally and retinal y superiorly in a right-eye chart. Positive tissue depth points from the vitreous toward the RPE. Visual y is inverted. The common chart uses the retained Watson A6 retinal-to-visual mapping and its numerical inverse; the independently fitted A5 and A6 transforms in the original engines remain unchanged for legacy reproduction. Device rotation and translation act on the physical coordinates before sampling the common scene or evaluating fields.

The displayed retinal surface has a 12 mm spherical radius and a smooth foveal depression. Layer thicknesses vary smoothly with eccentricity. At the center, the illustrative NFL/GCL/IPL/INL/OPL/ONL/IS/OS/RPE thicknesses are 8/3/3/5/15/70/24/30/12 µm. Away from the center they approach 20/42/46/39/15/45/24/30/12 µm. These values and the 340 µm central ganglion displacement are assumed display/model parameters, not an imported histological measurement set.

The microscopic view includes rods, cones, their inner/outer segments and nuclei, RPE mosaic, bipolar cells with ON/OFF IPL stratification, horizontal and amacrine processes, and ganglion somata, dendrites and axons. Macroscopic fibers reuse the existing published Jansonius reference. Cell density is reduced for inspection; 225, 529 or 961 population sites are computational samples, not estimates of anatomical cell count. Foveal photoreceptors and displaced inner cells are represented schematically.

Ganglion arbors are seeded, space-filling branching trees through a stratified target cloud, with a path-length penalty and tapered processes. They are procedural morphologies, not microscopy reconstructions. The selected cable uses those same physical nodes. The form is informed by retinal histology and the generative approach of Italiano et al.; no ModelDB reconstruction or source code is redistributed. Soma size, midget/parasol proportions and arbor diameter are assumptions. There is no individual patient anatomy.

References: Webvision, human fovea; Italiano et al. 2022, model and provenance; Watson 2014.

Three distinct physical pathways

PathwayIncident quantityPhysical calculationNeural input
Epiretinal electricalCommanded current, µABalanced finite disks in an infinite homogeneous isotropic conductorExtracellular potential differences, mV
Subretinal PRIMA880 nm irradiance, mW/mm²Two-series-diode photovoltaic circuit, electrode polarization, shared return and half-space fieldExtracellular potential differences, mV
Direct optogeneticWavelength and irradiance, mW/mm²Attenuation and opsin state kineticsLight-dependent conductance and photocurrent, pA

The shared scene is a luminance raster, sampled at the device’s actual location. Electrical encoding samples each contact aperture; optional edge encoding, channel gains and total-current budgeting are applied before field calculation. PRIMA projects the scene into its local 2 mm square and averages illumination over its reconstructed 378 hexagonal pixels. Pixel gains affect photocurrent, not the original projected image. The two recorded stages remain separate. Projection offsets are in the implant’s local frame. The native ProViSim image calculation remains conditional on projected input, not on the calculated photovoltaic current.

The retained electrical solver is full-space; PRIMA is a half-space approximation with its own shared-return interface equations. They are not interchangeable boundary conditions. Curvature is a registered display surface, not a curved conductive domain. Neither engine solves layered or anisotropic tissue, and the UI offers no misleading anisotropy control.

Electrical probes use soma and dendrite depths below the local NFL surface. PRIMA’s lower bipolar probe is at the local INL bottom, with the upper probe a configurable bipolar extent toward the vitreous. The device plane is one gap below that reference. The cable transfer uses each node’s depth relative to this same INL reference, all in the local planar approximation. This is an assumed subretinal registration, not a measured implant-to-cell reconstruction.

The electrical field overlay is a cathodic reference snapshot. Anodic-first selection reverses delivered waveforms and the active cable’s field; it does not introduce a validated polarity-dependent percept model. Interleaved trains are individually resolved in the detailed cable. The reduced population uses the selected-contact field snapshot and a common pulse envelope; it is not a full per-contact interleaved network calculation. All reported phase charges use complete biphasic pairs only.

PRIMA’s network input uses the final-pulse spatial field pattern multiplied by its calculated array-current envelope, normalized at the same snapshot. The detailed cable instead uses every captured pixel current and its shared-return transfer matrix. No optical irradiance is added to electrical current. Hybrid mode combines electrical fields and direct opsin conductance at the neural stage; optical delay allows sequential or overlapping delivery.

Response tiers and equations

Tier 1: bounded recruitment

The static surrogate is R(d,θ) = 1 − exp(−max(d,0)/θ). Ganglion and bipolar driving differences are divided by explicit threshold scales and weighted by survival. Opsin conductance and surviving natural input can contribute. The displayed result is peak unitless recruitment, without membrane voltages or spikes. This is a screening approximation, not an experimentally measured threshold distribution.

Tier 2: reduced retinal network

Each site has ON or OFF bipolar input, a low-pass bipolar state, lateral amacrine inhibition and a leaky integrate-and-fire ganglion state. Bipolar and amacrine states use exponential updates. The default bipolar time constant is 15 ms; lateral inhibition is 0.3 and the amacrine time constant 25 ms. Optional natural input is weighted by photoreceptor survival and an 80 ms adaptation state. ON sites rectify input minus half its adaptation; OFF sites rectify adaptation minus input. This is a simplified contrast circuit, not a full retinal microcircuit.

Ganglion voltage relaxes toward −65 mV + excitation with a 12 ms time constant. A −50 mV threshold resets voltage to −65 mV and imposes a 2 ms refractory interval. Excitation is surviving ganglion fraction times the rectified sum of electrical recruitment and direct/synaptic input. Default electrical and synaptic gains are 750 and 220 mV, respectively; both are exposed assumptions. They set a useful response range and are not fitted physiological measurements. Opsin conductance is 3 nS × expression × open fraction, with a 0 mV reversal, and its current is normalized by an assumed 180 pA scale for population coupling.

The response is spikes per second over the chosen observation window. Spikes, bipolar/amacrine states and voltages are exported. Threshold sensitivity at factors 0.5, 1 and 2 is a model sensitivity study, not a statistical confidence interval. Comparisons retain identical scene, tissue, sampling, seed and response parameters; device pulse schedules remain explicit.

Tier 3: selected branched active cable

The surrounding population remains Tier 2. One selected ganglion morphology is integrated with membrane capacitance, axial resistance, leak, Na, delayed K, A-type K, excitatory synaptic conductance and optional opsin conductance. A tree elimination solves implicit axial/membrane coupling; gates use Rush–Larsen updates. Axial resistance includes both adjacent half-sections. Extracellular coupling uses potential differences, so a constant extracellular offset cannot stimulate the cell.

Membrane capacitance is 1 µF/cm² and axial resistivity 110 Ω·cm. Maximum Na/K/A densities in S/cm² are dendrite 0.04/0.012/0.036, soma 0.07/0.018/0.054, initial segment 0.7/0.018/0.054, and axon 0.07/0.018/0. Reversal potentials are Na +35 mV and K −75 mV; leak is 0.000005 S/cm² at −62.5 mV. These follow the retinal Fohlmeister-family implementation documented in ModelDB 267391, with independently written integration. Calcium and calcium-activated K dynamics are omitted; this is not a reproduction of the complete published model. Temperature dependence and compartment-refinement convergence are not established.

There is 20 ms dark equilibration. Dendritic synaptic conductance follows the selected site’s calculated bipolar state using an assumed 0.002 nS/µm² scale. Direct opsin expression uses that scale per membrane area and the configured open fraction. The selected cable’s spike count is shown beside its reduced-network count to reveal disagreement. The population output decoder continues to use Tier 2; one detailed cell is not substituted for the entire retinal population.

References: Fohlmeister and Miller 1997; Italiano et al., ModelDB 267391.

Optogenetics and optical energy

Existing opsin modules are reused from the Optogenetics Laboratory: ChR2 uses its Williams-family kinetics; Chronos and ChrimsonR use explicitly reduced state models. Their spectral weights and kinetics are assumptions at the retinal expression site, not evidence of human retinal calibration. The clinical rationale for a distinct optogenetic pathway is described by Sahel et al. 2021.

Transmission to depth z is exp(−[μa + μs′(λ/590 nm)^−b] z), with depth converted to mm. This is a wavelength-scaled effective attenuation model. Reduced scattering is treated as removal from the modeled beam; there is no diffuse transport or lateral redistribution. A 4 mm square footprint, its placement and the common image determine which cells receive light. It is not Monte Carlo radiative transfer.

The thermal trace solves a lumped energy balance for a peak-illuminated tissue column: C_A dT/dt = I_abs − C_A T/τ, where C_A = volumetric capacity × thickness and I_abs = incident irradiance × [1−exp(−μa thickness)]. It uses assumed 0.004 J/mm³/K capacity, 250 µm thickness and 500 ms cooling by default. Scattering is not counted as absorbed heat. This is not a spatial temperature map, ocular exposure limit or safety determination. For PRIMA, illumination is at the implant; displayed heating is only a common illustrative absorption calculation, not heating predicted from the actual silicon package or retina.

Tissue state, vision and time

Disease presets alter photoreceptor, bipolar and ganglion survival, which changes the neural calculation as well as displayed cells. RP-like loss rises peripherally; geographic atrophy and focal injury use a smooth lesion region. Healthy tissue has unit survival. These are explanatory spatial masks, not histological fits or predictions of eligibility for an implant. Photoreceptor survival affects simulated activity only when surviving natural input is enabled. Direct ganglion stimulation can bypass photoreceptor loss; photovoltaic recruitment requires the assumed bipolar pathway.

The coupled percept is a fixed-scale Gaussian rendering of output samples. It is not acuity, reading performance, clinical vision or a cortical reconstruction. Electrical spatial and biphasic axon-map models and PRIMA’s ProViSim baseline are retained separately. The electrical baseline can show axonal elongation; the new Gaussian population decoder does not reproduce that axon-mediated percept. Baseline brightness units differ between technologies. PRIMA’s image-only baseline may be nonzero at zero irradiance, which is explicitly distinct from the zero coupled neural result.

Local images and videos are converted in the browser. Videos use up to eight center-cropped grayscale samples from the first two seconds. Moving inputs and video frames run as independent trials that restart the circuit and network. Playback within a trial shows the calculated response; sampled-sequence playback is not a continuous video-driven neural integration with state carried between frames.

Verification and remaining work

The original electrical and PRIMA test fixtures, upstream licenses, model versions and methods are retained. New tests cover shared image orientation, exact coordinate round trips, malformed imports, geometry, seeded anatomy, zero stimulus, lost cells, opsin expression, placement, charge balance and polarity; analytic optical transmission and lumped heating; cable SI units, dense-matrix residuals, passive RC dynamics, gauge invariance and time refinement; PRIMA field-transfer equivalence and passive current observation; and numerical regression of imported experiments. Graphics-enabled CI checks real workers, WebGL rendering, editing, comparison, legacy presets, exports and responsive layouts.

Biological validation against compatible retinal recordings and perceptual validation against matched patient outcomes remain unperformed. Also pending are authentic traced morphology imports, layered/anisotropic field solvers, calibrated diffuse optical transport, custom photovoltaic pixel geometry, continuous video integration, population-wide detailed electrophysiology and validated task-performance predictions. Circular custom electrical arrays are supported; unsupported Designer shapes or multiple planes are rejected rather than silently approximated.

Experiment JSON stores settings, seed, units, model versions, assumptions and calculated arrays. Import restores settings and recomputes; imported results are not trusted. The RTL- identifier is a non-cryptographic label, not a content-integrity or provenance guarantee. The original reference asset’s SHA-256 remains verified.

Retained electrical reference

The following original methods describe the unchanged electrical baseline and its original scope. Their statements about absent membrane or network models refer to that retained engine, not the exploratory unified tiers above.

Scientific scope

This laboratory connects editable electrode geometry and an exact pulse schedule to two distinct calculations: a finite-contact electrical field and a published phenomenological percept model. Position and delivered current are shared inputs. An electric-field magnitude is not converted into a phosphene through an invented activation threshold.

The implementation supports a pinned Argus-II reference geometry, custom regular and hexagonal arrays, both eyes, position, rotation, contact spacing and diameter, retinal gap, finite returns, disabled contacts, scene encoding, current budgets, and simultaneous or interleaved delivery. Seven examples isolate different assumptions. Pin one experiment before changing a parameter for a controlled comparison.

The percept model represents human epiretinal stimulation phenomenologically. The electrical model is an infinite homogeneous isotropic conductor. Retinal layers, anisotropy, implant insulation, electrode polarization, membrane voltage, spikes, synaptic processing, degeneration and learning are not solved. Outputs are not patient-specific predictions, visual acuity, recognition scores or electrochemical safety determinations. Subretinal photovoltaic and suprachoroidal devices require different models.

The 3D retina supplies anatomical context. Its curvature is a display assumption; electrical calculations use the explicitly labeled tangent plane. The fibers are population-average reference trajectories, not a measured individual retina.

Published percept models

The browser reproduces pulse2percept 0.11.0, including the circular spatial model, the Beyeler axon map, and the Granley–Beyeler representative biphasic model. The latter has a matched pulse-scaled circular comparison. Unmodified Python reference code generates the anatomy and independent golden outputs. The browser uses JavaScript and static local assets, with no remote scientific service.

For a candidate soma s, its soma-to-disc path A(s), electrode position e, and along-path distance dA(s,p), the spatial response is:

P(s) = max over p in A(s) of
       sum over electrodes e [ a_e × exp(−|p−e|²/(2ρ²) − dA(s,p)²/(2λ²)) ]

The electrode contributions must be summed at each path point before the maximum is selected. Adding individually calculated single-electrode percepts would produce a different model. The circular model uses Gaussian spatial spread without the axonal path-distance term. Beyeler et al., 2019.

The effective spatial spread ρ and axonal spread λ are in micrometers. They are adjustable percept parameters, which require percept measurements for patient-specific fitting. Electrode height, diameter and conductivity do not appear in these equations. Those physical inputs therefore change the electrical solution without silently changing fitted percept parameters.

Pulse-dependent representative image

For amplitude a in multiples of the supplied detection threshold, phase duration t in milliseconds and frequency f in hertz, the pinned package uses:

a_corrected = a × (0.054326 + 2.095t)
brightness factor = 0.1492147 × a_corrected + 0.0163851f
size factor       = max(1.0812 × a_corrected − 0.35338, 10²/ρ²)
streak factor     = max(1.56 − 0.54 × t^0.21, 10²/λ²)

The amplitude drive becomes the brightness factor, and the spatial and axonal variances become ρ² × size factor and λ² × streak factor. These coefficients reproduce the package’s documented Argus-II phase-duration threshold refit; they are not presented as the unchanged coefficients from the original paper. Granley & Beyeler, 2021.

The threshold anchor is the 50%-detection current at 0.45 ms phase duration and the selected train frequency. This interface supplies an assumed common threshold. Changing frequency does not create a new measured calibration. The brightness expression is not a detection probability, and no hard threshold gate has been added.

This is a representative pulse-train image, not a time-resolved movie. There is no modeled adaptation, interphase-gap effect or pulse-by-pulse dynamics. The pulse-dependent model requires simultaneous zero-delay trains and rejects interleaved delivery. A spatial model may compare interleaved electrode amplitudes only under an explicitly timing-independent assumption.

Anatomy, coordinates and geometry

The reference generates 1,000 Jansonius trajectories with 500 radial samples each. After unused paths are removed, 638 bundles and 123,122 path segments per eye remain. Both eyes retain their actual native coordinates and closest-path assignments. Display paths are thinned separately and are never solver input. The model uses minimum axonal sensitivity 1e-3, minimum current spread 1e-8, zero percept threshold and explicit meridian_blend=0. Jansonius et al., 2009.

Watson’s fitted coordinate transforms preserve x sign and invert y sign in the pinned package’s retinal basis. The physical plan has x increasing rightward and y upward; percept image row zero is the superior visual field at +15°. The separately fitted forward and reverse polynomials are not exact inverses. Both 0.5° and 0.25° visual grids span ±15°. Watson, 2014.

The editable Argus-II reference starts at 6 × 10 contacts, 575 µm pitch and 225 µm diameter, following pulse2percept 0.11.0. Published descriptions can use different dimensions. This label identifies a software reference geometry, not verification of a particular implanted device.

Electrical field

Each circular contact is a prescribed, uniformly injecting current sheet in an infinite conductor with conductivity σ. With signed current I_i, contact area A_i, and all lengths converted to meters:

φ(r) = sum_i [ I_i/(4πσA_i) × integral over disk_i (1/|r−r'|) dA' ]
E(r) = −gradient φ(r)
J(r) = σ E(r)

Positive current injects into the conductor. A cathodic stimulating contact has negative current; finite return disks have positive current. All source and return contributions are included, even outside the plotted region. Potential is relative to infinity. Outputs are extracellular potential in mV, electric field in V/m, and bulk current density in A/m².

For a single disk of radius a, axial separation h > 0 and current I, the exact full-space potential provides an independent closed-form check:

φ(h) = I / [2πσ × (sqrt(a²+h²) + h)]

A +80 µA source with 100 µm radius, 100 µm separation and 0.2 S/m conductivity gives +263.69654 mV and −1864.61614 V/m axial field below the disk. A finite return changes the total; these values are the isolated-source contribution.

The source sheet has conductor on both sides. It is not an equipotential metal contact or an electrode on an insulating half-space boundary. Those boundary conditions produce different solutions. The default 0.2 S/m conductivity is an editable demonstration assumption. A layered anisotropic retina and an electrode interface need separate parameters and verification. Esler et al., 2018; Chen et al., 2020.

Numerical method

The default disk integration uses 32 Gauss–Legendre radial nodes and 128 periodic angular nodes. Small separation/radius ratios automatically raise the order. Source-plane observations and ratios below 0.025 are rejected. Cubic Hermite interpolation preserves the relationship between radial field and potential derivative. Geometry-dependent transfer coefficients use a bounded 128 MiB cache; large requests recompute evicted coefficients without coarsening the grid.

Independent Python references integrate azimuth analytically with complete elliptic integrals, then adaptively integrate radius. This differs from the JavaScript quadrature. Closed-sphere current-flux, gradient, polarity, conductivity-scaling and finite-return tests supplement pointwise comparisons.

Delivery

The scheduler stores exact piecewise-constant cathodic and anodic intervals, including the interphase gap. It delivers only complete biphasic pulse pairs inside the selected train window. Charge is not estimated by rounding edges to time bins.

Phase charge [nC] = |current| [µA] × phase duration [ms]
Geometric phase density [µC/cm²] = 100000 × phase charge [nC] / (πa² [µm²])

For example, 100 µA for 0.2 ms delivers 20 nC per phase. A 100 µm-radius contact has geometric phase density 63.66198 µC/cm². The ideal symmetric waveform has zero net charge, within floating-point arithmetic. Total absolute charge counts both phases of delivered stimulating-contact pulses, excluding return charge. Arithmetic charge balance alone does not establish electrochemical safety.

Interleaved contacts use equally spaced slots at the selected per-contact repetition frequency. The scheduler rejects pulses that cannot fit. The electrical map is a named cathodic snapshot: all simultaneous contacts, or one interleaved contact. Contact-site probes use the same snapshot.

Encoding and current budgets

Synthetic letter, Landolt C, bar and checkerboard inputs are defined in visual degrees. Intensity encoding averages a 3 × 3 aperture; edge encoding uses a clipped Sobel response. The aperture half-step is max(0.15°, pitch/280/4). The 280 µm/degree factor is an explicit encoder convention, not the Watson coordinate transform. Single- and two-contact probes bypass the scene encoder.

Unavailable contacts are zeroed before current-budget scaling. Scene weights multiply a per-contact peak-current limit. The optional summed-current budget scales all amplitudes down together when their sum exceeds the cap. It does not increase weak channels to force expenditure of the cap. Summed channel current differs from instantaneous current in an interleaved slot.

Measurements

OutputDefinition
Percept brightnessRaw model units, without independent normalization of comparison images. Different model families are not automatically luminance-calibrated.
Half-peak support areaPixel count at ≥50% of image peak times sample area, in deg²; not an exact contour integral.
Support centroid and elongationUnweighted support geometry; moment-axis widths are four times square roots of covariance eigenvalues.
Integrated brightnessSum of image samples times visual sample area; depends on grid and field of view.
ComponentsFour-connected support-pixel regions, not a count of separately perceived phosphenes.
Electrical extremaSampled grid maxima, not certified continuous-space maxima.
Contact probesField and potential evaluated at the actual contact-center coordinates on the retinal plane, separately from raster sampling.
Zero-drive resultsExact zero images and fields, with undefined centroids and elongations represented as unavailable.

Electrical maps offer 61, 121 or 241 samples per side over the array or a local contact-centered region. A wide grid can miss a narrow peak between samples. Boundary-tail flags and sampling warnings expose some limitations; absence of a warning is not a convergence certificate.

Changing an input clears stale computed results. A pinned trial retains its settings and numerical output. Versioned JSON stores settings, full samples, waveforms, contact probes, model identity and reference checksum. Restoration validates the definition and recomputes rather than trusting imported results. CSV headers carry coordinates and units; initial # lines contain provenance.

Verification

Fresh recovered-build checks include 16 independent native percept cases, 69 independent field fixtures, 29 percept tests, 9 JavaScript field tests, 15 complete-experiment/client tests, and 9 Python field tests. The 62 targeted scientific tests are distinct from the repository build and browser checks. All 52 generated percept assets and fixtures reproduced byte-for-byte in a separate output directory.

Numerical checkObserved recovered-build result
Percept implementation maximum error / reference peak2.242 × 10⁻⁷
Percept implementation RMS error / reference peak2.743 × 10⁻⁸
Direct field relative potential error4.631 × 10⁻⁸
Direct field relative vector-E error5.414 × 10⁻⁶
Cached field relative potential error4.472 × 10⁻⁸
Cached field relative vector-E error9.650 × 10⁻⁶
Visual-grid refinement integral change, three cases0.0285–0.6406%
Doubled trajectory resolution, peak-normalized RMS change0.5578–0.9648%
Doubled trajectory resolution, largest local change / peak18.09–26.84%

These are measured errors on specified cases, not universal bounds. The native trajectory study still shows substantial local differences near nearest-bundle transitions. Uniform pointwise convergence and patient validation have not been established. Close implementation agreement does not remove uncertainty in the underlying model or its fit to an individual.

The largest electrical stress case retains all 58,081 points on a 241² grid and stays within the cache limit. When the transfer set exceeds the cache, repeat runs can remain slower than smaller arrays. Execution times are host-dependent; the worker remains cancellable.

Reproduce

npm ci
node --test tests/retinal-field.mjs tests/retinal-percept.mjs tests/retinal-experiment.mjs
npm run build --ignore-scripts
node scripts/verify-retinal-build.mjs
python3 -m venv /tmp/retinal-reference
/tmp/retinal-reference/bin/pip install -r scripts/retinal-reference/requirements.txt
/tmp/retinal-reference/bin/python scripts/retinal-reference/generate.py
/tmp/retinal-reference/bin/python -W error scripts/retinal-reference/field_reference.py
/tmp/retinal-reference/bin/python -m unittest discover -s scripts/retinal-reference -p 'test_field_reference.py'
/tmp/retinal-reference/bin/python scripts/retinal-reference/convergence.py

The anatomy is stored as 31 small compressed chunks. Lossless corresponding-eye XOR, byte-plane rearrangement and gzip reduce transfer size; the loader checks every chunk, reverses the transforms and verifies the complete native-data SHA-256. No scientific coordinates or arc distances are rounded by compression. The full upstream BSD-3-Clause license is included.

Sources and next scientific steps

The primary papers linked above specify the model families and their context. The pinned pulse2percept release and committed generators define the numerical reference. Official model and coordinate documentation provide additional background, but live documentation may change.

A physical extension should specify tissue layers, tensor conductivity, implant boundaries and interface dynamics before adding neural activation. A neural extension requires a separately verified membrane model and consistent extracellular coupling. Patient-specific evaluation needs measured anatomy, electrode position, thresholds and percepts, with parameter fitting separated from held-out evaluation. See Kish, Yuan & Weiland, 2023 for the additional anatomy and neural modeling involved in that goal.