The advent of large high-dimensional datasets in neuroscience has been an important milestone for advancing our understanding of neural information processing and improving performance of brain computer interfaces. However, most existing methods of analysis fall short of capturing the complexity of concerted population activity. Novel techniques need to address this complexity and be applicable in a wide range of neural data analysis scenarios. In this work, we built non-parametric pair copula constructions which disentangle single-neuron statistics from the dependency structures within the population. This approach makes it possible to study dependencies between variables with vastly different statistics (e.g. continuous behavioural variables and discrete spike-trains). We estimate copula and margin densities with Normalizing Flows (NFs) to avoid misspecification of parametric models, especially in the case of discrete variables. Overall, NFs achieved comparable performance to existing non-parametric estimators when trained on artificial data with known dependency structures, while allowing for easier sampling and more flexibility. Finally, we show our framework’s aptitude to capture non-symmetric dependencies in groups of neurons in the rodent primary visual cortex responding to a visual learning task.