Neuronal network reconstruction is an important problem in neuroscience as it helps understanding neuronal circuit and function. With the advancements of the calcium imaging technique, the dynamic activity of hundreds of neurons can be observed, which also provides a foundation for modelling and inferring network connectivity directly from data. In this paper, the dynamic behavior of each neuron is described using a non-linear autoregressive exogenous input (NARX) model. Since the neurons are inter-connected, the NARX models contain inputs which depend on outputs of the surrounding neurons. The complete model is expressed in a linear-in-parameter form and the network dynamics are identified by obtaining a block-sparse solution given the calcium oscillation data. The solution to the identification problem also provides the network topology simultaneously, i.e., the physical connections between the neurons. To demonstrate the accuracy of the proposed method, two experimental case studies are considered and the results are compared with other network identification methods. © 2019, IFAC (International Federation of Automatic Control) Hosting by Elsevier Ltd. All rights reserved.