Multi-electrode arrays (MEAs) can record extracellular action potentials (also known as 'spikes') from hundreds or thousands of neurons simultaneously. Inference of a functional network from a spike train is a fundamental and formidable computational task in neuroscience. With the advancement of MEA technology, it has become increasingly crucial to develop statistical tools for analyzing multiple neuronal activity as a network. In this paper, we propose a scalable Bayesian framework for inference of functional networks from MEA data. Our framework makes use of the hierarchical structure of networks of neurons. We split the large scale recordings into smaller local networks for network inference, which not only eases the computational burden from Bayesian sampling but also provides useful insights on regional connections in organoids and brains. We speed up the expensive Bayesian sampling process by using parallel computing. Experiments on both synthetic datasets and large-scale real-world MEA recordings show the effectiveness and efficiency of the scalable Bayesian framework. Inference of networks from controlled experiments exposing neural cultures to cadmium presents distinguishable results and further confirms the utility of our framework.