Simulation Based Inference of a Simple Neural Network Structure

arXiv:2604.18599 · stat.AP, math.ST, q-bio.NC, stat.TH · Submitted 2026-04-09 · Read on arXiv

Listen

Radio episode about this paper

Transcript

Introduction to the show: ident: Genomics Radio. Generated commentary on the latest computational biology and genomics papers.

Ines: Today's paper: "Simulation Based Inference of a Simple Neural Network Structure".

Marcus: The gist: The proposed method uses simulation-based inference to estimate parameters of an underlying random graph model from sparse spike train statistics, offering better results than traditional sub-network reconstruction methods.

Ines: First, who's behind it and why it matters.

Title and authors: Ines: So this paper is called Simulation Based Inference of a Simple Neural Network Structure. I mean, it sounds very technical, but basically, they are trying to figure out the structure of a whole neural network when they only have data from just a few recorded neurons.

Marcus: Right. They're dealing with that huge problem where you can only see a tiny fraction of the actual neurons in the brain. This paper proposes using simulation to estimate those parameters of the underlying random graph model without having to reconstruct the whole thing directly, which is what traditional sub-network reconstruction methods try to do.

Yuki: From my side, I’m thinking about how this relates to population genetics. If you only sample a small group of individuals, can you really tell anything about the entire species' evolutionary structure? This method seems to be looking at that scaling issue in a different way by focusing on simple spike train statistics like empirical frequency or interspike intervals.

Ines: Exactly. What they are doing is taking their numerical model of the whole network, which has parameters like the number of neurons and connection probability, and then generating artificial data from it to mimic real experiments. They sample that artificial data in a way that reflects how sparse the real data is supposed to be.

Marcus: So instead of trying to solve for whatever is most likely given the exact observed spike trains, they are solving for the parameter p of a random graph model by maximizing something else based on those simulated statistics. It's about finding that connection probability p, which governs how connections are made in this network.

Yuki: And it seems like they're aiming to address that massive under-sampling problem by reproducing it numerically as faithfully as possible, which is what they call the essence of their simulation based inference approach <ref:2604.18599#pg3>.

Ines: That’s the core idea, isn't it? They are using a numerical model of the whole network to generate artificial data and then sampling that in a way mimicking real experiments to adjust their model parameter p <ref:2604.18599#pg3>.

Marcus: And they change how they look at the problem itself. Instead of solving for something like arg max p in P Y = y which is what some other methods do, they are aiming for arg max p in Pp T(Y) = T(y) <ref:2604.18599#pg10>.

Ines: And the way they define that statistic T is key here. It allows them to use spike train statistics, like the empirical spikes frequency or the interspike interval distribution, to estimate p <ref:2604.18599#pg1> and see if those later depend critically on the random graph parameters <ref:2604.18599#pg3>.

Marcus: It’s interesting that they are moving away from models that ignore anatomical knowledge, which they say is a major barrier to dialogue between experimentalists and statisticians <ref:2604.18599#pg3>.

Yuki: And for someone looking at the big picture, it suggests that by focusing on these simpler statistics, we might be able to make inferences about the larger network structure even when the observed sample is tiny.

The paper's summary: Ines: So what’s actually recovered from this approach? The paper summarizes it by saying they focus on simple spike train statistics like the empirical spikes frequency or the interspike interval distribution, which are much easier to estimate than trying to reconstruct the whole network structure directly.

Marcus: They argue that these simpler statistics can be estimated by simulations, and if you have a few observed spike train statistics, they provide enough information to estimate p <ref:2604.18599#pg1>. This is the main payoff for using this simulation based inference approach <ref:2604.18599#pg3>.

Yuki: From a population genetics viewpoint, if we only see a small subset of observed spike trains, can these simple statistics really support conclusions about the underlying graph parameters? I mean, what’s the actual power of this inference on that level?

Ines: They are testing how sensitive those later statistics are to the random graph parameters <ref:2604.18599#pg3>. They want to see if they actually depend critically on p <ref:2604.18599#pg3>. This is their central working hypothesis, which is that the graph is a realization of a directed Erd˝os-R´enyi model where P(Wi→j = w) equals p and P(Wi→j = zero) equals one minus p <ref:2604.18599#pg4>.

Marcus: They set up the model with parameters like the number of neurons and the connection probability between any pair of neurons, which are the basic ingredients they use <ref:2604.18599#pg3>. They also have random synaptic weights Wj→i in zero or one, where P(Wi→j = w) equals p and P(Wi→j = zero) equals one minus p with independent weights <ref:2604.18599#pg4>.

Yuki: So they are modeling the whole system—the neurons, their interactions, the spike generation—all numerically. That sounds like a lot of work to set up just to see if it actually yields better results than looking at just the recorded neurons <ref:2604.18599#pg2>.

Ines: They describe the dynamics using a discrete time model where Xt(i) equals one if there’s a spike in window t for neuron i, and Vt(i) is defined by summing contributions from all presynaptic neurons j based on their last spiking time Lt(i) <ref:2604.18599#pg5>.

Marcus: The core of the dynamics is that a neuron decides to spike at time t plus one with probability phi i(Vt(i)), independently of the others <ref:2604.18599#pg8>. This translates into a recursive way where Xt+one(i) equals one if Ut+one(i) less than or equal to phi i(Vt

i: ) and zero otherwise <ref:2604.18599#pg4>.

Yuki: When you boil it down, they are taking these spike train dynamics, which are really complex temporal processes, and trying to connect them back to the underlying network structure using these simpler statistics <ref:2604.18599#pg10>. It’s a big step in bridging the gap between observation and theory.

Ines: Right. They are not just looking at what the recorded neurons do, but they are using this simulation to infer the structure that generated those observations, which is what makes it different from traditional methods <ref:2604.18599#pg2>.

The paper's improvements: Marcus: So where does this approach improve things? The authors suggest they are addressing the massive under-sampling problem by reproducing it numerically as faithfully as possible, which is how they estimate the sampling distribution of their statistics <ref:2604.18599#pg3>.

Ines: They also propose a way to transform the original maximization problem. They shift it from solving arg max p in P Y = y to solving arg max p in Pp Y s/j=one Pp phi(Yj) = phi(yj) <ref:2604.18599#pg8>.

Marcus: That change means they can use the empirical spike frequency as an example of the function phi, which is just a count of spikes per time step <ref:2604.18599#pg10>. They generate K equals fifty independent directed graphs for each p to estimate Pp phi(Y) = phi(y) <ref:2604.18599#pg12>.

Yuki: Generating fifty graphs for every value of p sounds computationally intensive, but the point is they are trying to get a robust estimate of the probability distribution related to those statistics <ref:2604.18599#pg12>.

Ines: And once they have that estimate, they define their final estimate of p as the location of the maximum of a quadratic polynomial interpolating Qs/j=one Pbp T(Yj) = T(yj) between two points around the estimated value <ref:2604.18599#pg28>.

Marcus: That interpolation step is what helps overcome some bias, using Lemma two which gives an efficient algorithm to compute this maximizing value <ref:2604.18599#pg28>. It sounds like they are refining the estimate of p by looking at local behavior in the sampling distribution <ref:two thousand six hundred four point one eight five nine nine#pg28.

Yuki: And they also suggest that using the inter-spike interval, which is a compression of a neuron’s spike train, can be a useful statistic because its distribution can be seen as a vector in R100 on which dimension reduction via PCA yields a curve well parameterized by p <ref:2604.18599#pg25>.

Ines: So it’s like they are taking this complex spike train data, summarizing it down into something simpler that has a clear relationship with p <ref:2604.18599#pg27>. They propose using the estimated alpha parameter from fitting the ISI distribution to a Gamma distribution as a good statistic to distinguish between different values of p <ref:2604.18599#pg26>.

Marcus: And they also show that even when they remove reciprocal connections, they can still estimate the relevant connection probability p prime, which is defined as p prime equals p minus p squared over two <ref:2604.18599#pg13>. That’s a good sign, because it shows that their fitted model doesn't necessarily need to exclude reciprocal connections <ref:2604.18599#pg19>.

Yuki: That implies that the underlying structure they are finding is robust enough even with these kinds of assumptions about how the network is wired.

Conclusion: Ines: So to wrap up, this paper, "Simulation Based Inference of a Simple Neural Network Structure," proposes using simulation-based inference to estimate parameters of a simple neural network structure from sparse spike train statistics. The main implication is that they can use these statistics—like empirical frequency or interspike interval distribution—to get an estimate for the connection probability p <ref:2604.18599#pg1>.

Marcus: They show that their method, despite its complexity, gives results similar to classical graph reconstruction methods, though the optimal mean absolute error and standard error are about twenty times larger than those of the classical graph reconstruction method <ref:2604.18599#pg19>.

Yuki: But they also provide confidence intervals for p by assuming a distribution D equals minus two log Pbp0 T(Y) = T(y) minus log Pbpb T(Y) = T(y) is approximately chi squared one <ref:2604.18599#pg28>. So, they can actually put a measure of uncertainty around their estimate of p <ref:2604.18599#pg19>.

Ines: And the heuristic approach with the inter-spike interval distribution, using its shape parameter alpha derived from fitting it to a Gamma distribution is shown to be a good statistic for distinguishing between different values of p <ref:2604.18599#pg26>.

Marcus: Overall, they built this simulation based inference approach to tackle the under-sampling problem by numerically reproducing it as faithfully as possible <ref:2604.18599#pg4>. It’s a method that aims to provide better results than traditional sub-network reconstruction methods <ref:2604.18599#pg19>.

Yuki: It suggests that even with limited data, we can make inferences about the larger network structure by focusing on these simple statistics, which is something relevant for understanding evolution across species <ref:2604.18599#pg27>.

Ines: So that’s what they are doing with "Simulation Based Inference of a Simple Neural Network Structure." It’s about using numerical modeling to bridge the gap between sparse observations and the underlying network parameters.

Marcus: A lot of work on simulation-based inference is going into making these methods more efficient, and that's where we’ll be looking next.

Pierre Charitat, S. Geffray, C. Pouzat

IRMA University of Strasbourg and CNRS UMR 7501

stat.AP, math.ST, q-bio.NC, stat.TH

Submitted: 2026-04-09

Updated: 2026-10-14

Journal ref: Stochastic Processes and their Applications, 2026, 203 (105112)

License: http://arxiv.org/licenses/nonexclusive-distrib/1.0/

Importance score: 73/100

The gist: The gist: The proposed method uses simulation-based inference to estimate parameters of an underlying random graph model from sparse spike train statistics, offering better results than traditional

Key concepts

Simulation Based Inference (SBI)
This involves creating a numerical model of the entire neural network, generating artificial data from it under different parameters (like connection probability 'p'), and then sampling that artificial data in a way that mimics real experiments. The goal is to find the parameter 'p' that makes the simulated data match the observed spike train statistics.
Directed Erdős-Rényi Graph
This is a mathematical model representing a network where neurons are nodes and connections are directed edges. In this specific model, every possible connection between two neurons has an independent probability 'p' of existing, meaning it's a random graph structure.
Spike Train Statistics (e.g., ISI Distribution)
Instead of looking at individual spikes, the method uses summary statistics from the entire spike train, such as the empirical frequency or interspike interval (ISI) distribution. These statistics are used because they are more robust than raw data and can be effectively used to infer the underlying network structure.

Summary

The gist: The proposed method uses simulation-based inference to estimate parameters of an underlying random graph model from sparse spike train statistics, offering better results than traditional sub-network reconstruction methods.

Simulation Based Inference Approach

The core idea is to use a numerical model of the whole network to generate artificial data, which is then sampled in a way mimicking real experiments to adjust the model parameter p<ref:2604.18599#pg3>. This approach addresses the massive under-sampling problem by reproducing it numerically as faithfully as possible<ref:2604.18599#pg4>. The method involves solving arg max p∈(0,1) Pp T(Y) = T(y) instead of arg max p∈(0,1) Pp Y = y<ref:2604.18599#pg10>.

Data Generation Model

The model considers a directed Erd˝os-R´enyi graph involving two parameters: the number of neurons and the connection probability between any pair of neurons/vertices<ref:2604.18599#pg3>. The basic ingredients include a finite set I of n neurons, a set of random synaptic weights Wj→i ∈ 0, w, and a family of spiking probability functions ϕi<ref:2604.18599#pg2>. The central working hypothesis is that this graph is a realization of a directed Erd˝os-R´enyi model where P⁡(Wi→j = w) = p and P⁡(Wi→j = 0) = 1 − p, with independent weights<ref:2604.18599#pg4>.

Spike Train Dynamics

The discrete time model describes the spiking activity of a finite set I of neurons over time, where Xt(i) = 1 indicates the presence of a spike within the time window of index t<ref:2604.18599#pg5>. The membrane potential Vt(i) is defined by summing contributions from all presynaptic neurons j ∈ V·→i based on their last spiking time Lt(i)<ref:2604.18599#pg7>. A neuron decides to spike at time t+1 with probability ϕi(Vt(i)), independently of the others<ref:2604.18599#pg8>. The dynamics are formally translated into a recursive way where Xt+1(i) = 1 if Ut+1(i) ≤ ϕi(Vt(i)) and 0 otherwise<ref:2604.18599#pg4>.

Summary Statistics and Inference

The problem is transformed from solving arg max p∈(0,1) Pp Y = y to solving arg max p∈(0,1) Ys/j=1 Pp φ(Yj) = φ(yj)<ref:2604.18599#pg8>. A first example for the function φ is the empirical frequency of spikes per time step<ref:2604.18599#pg10>. The method involves discretizing p into a set P and generating K=50 independent directed graphs for each p to estimate Pp φ(Y) = φ(y)<ref:2604.18599#pg12>. The estimate ˆp of p is then defined as the location of the maximum of the quadratic polynomial interpolating Qs/j=1 Pbp T(Yj) = T(yj) between ˜p and its two nearest neighbors<ref:2604.18599#pg28>.

Empirical Results

The experiments show that variants of our method give very similar results, except for the alphas gaussian likelihood method, for small values of p<ref:2604.18599#pg17>. The optimal mean absolute error and standard error are about 20 times larger than the classical graph reconstruction method<ref:2604.18599#pg19>. The method can also provide confidence interval for the value of p by assuming D = −2 log Pbp0 T(Y) = T(y) − log Pbpb T(Y) = T(y) ≈ χ2 (1)<ref:2604.18599#pg28>. The empirical result shows that the estimated α parameter is indeed a good statistic to distinguish between the values of p<ref:2604.18599#pg27>.

Heuristic Approach to Find Statistics

The paper proposes using the inter-spike interval (ISI) distribution of a neuron as a relevant statistic because it is a compression of a neuron’s spike train<ref:2604.18599#pg24>. An ISI distribution can be seen as a vector in R100 on which dimension reduction via PCA yields a curve well parameterized by p<ref:2604.18599#pg25>. The estimated α parameter, derived from fitting the ISI distribution to a Gamma distribution, is shown to be indeed a good statistic to distinguish between the values of p<ref:2604.18599#pg26>. The average interspike interval length remains approximately constant for a given p, suggesting that using the estimated mean of the ISI distribution as the statistic on the spike trains is effective<ref:2604.18599#pg27>.

Removing Reciprocal Connections

The method was run with a neural network having reciprocal connections removed, and it still correctly estimates the relevant connection probability p′, which is defined as p′ = p − p2/2<ref:2604.18599#pg13>. This demonstrates that the fitted model does not necessarily need to exclude reciprocal connections<ref:2604.18599#pg19>.

Confidence Intervals

The method can also provide confidence interval for the value of p by assuming D = −2 log Pbp0 T(Y) = T(y) − log Pbpb T(Y) = T(y) ≈ χ2 (1)<ref:2604.18599#pg28>. This leads to an approximated 100α% confidence interval<ref:2604.18599#pg19>.

Numerical Details

The implementation of SBI methods requires efficient codes for simulations, and the authors developed their codes from scratch in Fortran<ref:2604.18599#pg29>. They used the xoshiro256++ generator of [3] to ensure reproducible results regardless of the compiler used<ref:2604.18599#pg29>. All our codes and scripts carrying out the simulations, their analysis and the generation of the figures of this article are available from our dedicated GitLab repository<ref:2604.18599#pg29>.

Acknowledgments

This article is dedicated to the memory of Antonio Galves, a brilliant mathematician whose curiosity could not be sated, and a very dear friend of the author (CP) who had the privilege to know him well<ref:2604.18599#pg31>. Pierre Charitat was supported by a fellowship from Strasbourg University Interdisciplinary Thematic Institute: Research in Mathematics, Interactions and Applications (IRMIA++)<ref:2604.18599#pg31>. S´egolen Geffray and Christophe Pouzat were supported by an ANR grant: SIMBADNESTICOST ANR-22-CE45-0027<ref:2604.18599#pg31>.

References

[1] Peter Barth´o, Hajime Hirase, Lena¨ıc Monconduit, Michael Zugaro, Kenneth D. Harris, and Gy¨orgy Buzs´aki. Characterization of neocortical principal cells and interneurons by network interactions and extracellular features. Journal of Neurophysiology, 92(1):600–608, July 2004<ref:2604.18599#pg32>.

[2] B. R. Bhat. On the method of maximum-likelihood for dependent observations. Journal of the Royal Statistical Society. Series B (Methodological), 36(1):48–53, 1974<ref:2604.18599#pg32>.

[3] David Blackman and Sebastiano Vigna. Scrambled linear pseudorandom number generators. ACM Transactions on Mathematical Software, 47(4):1–32, September 2021<ref:2604.18599#pg32>.

[4] Valentino Braitenberg and Almut Sch¨uz. Cortex: Statistics and Geometry of Neuronal Connectivity. Springer Berlin Heidelberg, 1998<ref:2604.18599#pg32>.

[5] D. R. Brillinger. Maximum likelihood analysis of spike trains of interacting nerve cells. Biological Cybernetics, 59(3):189–200, August 1988<ref:2604.18599#pg32>.

[6] David R. Brillinger, Hugh L. Bryant, and Jos´e P. Segundo. Identification of synaptic interactions. Biol. Cybern., 22(4):213–228, December 1976<ref:2604.18599#pg32>.

[7] George Casella and Roger L Berger. Statistical Inference. Second Edition. DUXBURY ADVANCED SERIES, Duxbury, 2002<ref:2604.18599#pg32>.

[8] E. S. Chornoboy, L. P Schramm, and A F Karr. Maximum likelihood identification of neural point process systems. Biological Cybernetics, 59(4–5):265–275, September 1988<ref:2604.18599#pg32>.

[9] Kyle Cranmer, Johann Brehmer, and Gilles Louppe. The frontier of simulation-based inference. Proceedings of the National Academy of Sciences, 117(48):30055–30062, May 2020<ref:2604.18599#pg32>.

[10] E. De Santis, A Galves, G Nappo, and M Piccioni. Estimating the interaction graph of stochastic neuronal dynamics by observing only pairs of neurons. Stochastic Processes and their Applications, 149:224–247, July 2022<ref:2604.18599#pg32>.

[11] Peter J Diggle and Richard J Gratton. Monte carlo methods of inference for implicit statistical models.

Improvements for AI systems

  1. Focus on inferring underlying network structure rather than reconstructing sub-networks by utilizing simple spike train statistics, like the empirical spikes frequency, or the interspike interval distribution, which can be estimated by simulations to provide information on the connection probability of the original network.

  2. Implement a simulation-based inference approach where one uses a numerical model of the whole network to generate artificial data and then samples spike trains, mimicking real experiments, allowing for an estimate of the parameters of a postulated random graph model.

  3. Use inter-spike interval (ISI) distribution as a summary statistic because our heuristic approach to find statistics on spike trains is to propose to use the inter-spike interval (ISI) distribution of a neuron, noting that it can be summarized by estimating its shape parameter, denoted as α of a fitted Gamma distribution.

  4. Develop dimension reduction techniques, such as Principal Component Analysis (PCA), on the ISI distributions generated for different values of the connection probability to project them onto a one-dimensional space to distinguish the associated value of p.

  5. Employ quadratic polynomial interpolation between discrete points in the sampling distribution of the statistic to overcome bias, using Lemma 2 gives an efficient algorithm to compute this maximizing value as described in Section C.2.

Related papers