Application of Bagged Copula-GP: Confirming Neural Dependency on Pupil Dilation
Authors: Maximilian Walden
Venue: Accepted by TMLR ()

Abstract: Advancements in recording techniques has enabled the ability to record thousands of neurons simultaneously, shifting the needs within the field of computational neuroscience to powerful computational and statistical techniques. Copula-GP is a recently developed state-of-the-art parametric mutual information estimator found to outperform other novel non-parametric methods when utilized on highly dimensional data. Here, we utilized Copula-GP together with Gaussian Process Factor Analysis (GPFA) to investigate the information interaction between neuronal processes within the visual cortex of live mice and pupil dilation. We found usage of GPFA as a preprocessing step to Copula-GP was an effective means of investigating neuronal dependence, allowing flexibility in analysis and finding results in agreement with prior literature, and additionally extended Copula-GP with a bagging framework, allowing for the aggregation of model estimations and allowing for more accurate estimation accuracy and representation of dependency shape. We validated our bagging algorithm on simulated data sampled from known distributions, and utilized bagged Copula-GP with GPFA on said neuronal data to find results in agreement with baseline Copula-GP but with more stability.

Published in Transactions on Machine Learning Research (12/2024)
Application of Bagged Copula-GP: Confirming Neural De- pendency on Pupil Dilation
Maximilian Walden University of Edinburgh
M. Walden @sms. ed.
Reviewed on OpenReview: https: openreview. net/ formm? id= hQIs1CHFGN
Abstract
Advancements in recording techniques have enabled the ability to record thousands of neu- rons simultaneously, shifting the needs within the field of computational neuroscience to that of powerful computational and statistical techniques. Copula-GP is a recently devel- oped state-of-the-art parametric mutual information estimator found to outperform other novel non-paramelric methods when utilized on highly dimensional data. Here, we utilized Copula-GP together with Gaussian Process Factor Analysis (GPFA) to investigate the in- formation interaction between neuronal processes within the visual cortex of live mice and pupil dilation. We found usage of GPFA a preprocessing step to Copula-GP was an effective means of investigating neuronal dependence, allowing flexibility in analysis and finding results in agreement with prior literature. We additionally extended Copula-GP with a bagging framework, allowing for the aggregation of model estimations and allow- ing for more accurate estimation accuracy and representation of dependency shape. We validated our bagging algorithm on simulated data sampled from known distributions, and utilized bagged Copula-GP with GPEA on said neuronal data to find results agreement with baseline Copula-GP GP but with more stability.
1 Introduction
Computational Neuroscience aims to answer how coguition, behavior, and learning are encoded in the brain through the use of machine learning and statistical methods on large amounts of neuronal data (Blundell ett al. 2018 McDougal et al., 2017 Vu et al. 2018), with answers to these questions having applications in fields such as prosthetics, medicine, and machine lcarning (Jangid et al. 2021). Some of the primary questions on this frontier include: How is information encoded in populations neurons, and how this information related behaviors observed (Kudryashova al.] (022)?
As advances in recording techniques have allowed for the recording of hundreds to thousands of neurons 2017 Dombeck et al.] 2007), the bottlenecks keeping researchers from answering these central questions in neuroscience have shifted away from data-related issues to the purely computational and statistical (Vu et al. 2018). Answering such questions with machine learning comes down to analysis of the intricate dighly-dimensional multivariate dependencies (i.e. dependencies and correlations between many random variables, many of which have different distrributions) in recorded neuronal and behavioral data, much of which varies across spatial and temporal dimensions (Kudryashova et al.) 2022 Jun al., 2017 Glaser et al. 2019; Quian Quiroga & Panzeri 2009). This becomes especially challenging
at once (Jun et al.
one considers
how different behaviors occur occur on scale of milliscconds (Pakan et al. 2018 Dombock et al.) 2007 Mathis et al. 2018), even more so as more variables are included in the data and the dependencies become even more computationally intense to analyze accurately (the "curse of dimensionality") (Koppen) 2000). Thus, in order to properly analyze such copious amounts of intricate high-dimensional data, we need methods that are both resistant to high dimensionality as well as equipped with the ability to properly analyze dependencies between temporally- disparate variables in time series.
scale of hours or days, whereas related neuronal activity

---

Published in Transactions on Machine Learning Research (12/2024)
Porametric Copula- GP model for analyzing multidimensional neuronal and behavioral relationships (Kudryashova et al. 2022) aimed to fulfill this requirement via novel framework for parametric estimation of multivariate distributions and mutual information through the use of copulas, probabilistic structure adept at mutual information estimation. The core paper's contributions to the body of work includes the development of the Copula-GP python package, a package which allows for the training of copulas, model selection, and copula. entropy extraction. Copula-GP was found to outperform other non-parametric mutual information estimators when predicting mutual information for a larger number of variables (2 10). This combined with it's robustness to the "curse of dimensionality" makes it well suited for analysis of highly dimensional neuronal data.
Here, we validated usage of Copula-GP together with Gaussian Process Factor Analysis (GPFA; a dimen- sionality reduction method robust to qualities of neuronal data through which we are able to extract latent variables), doing so through confirming statistical dependence of neuronal data extracted from mice on pupil dilation. Within the literature, there is very little use of GPFA with mutual information estimators, much less conditional parametric copulas on neuronal data (or on any data). This project aimed to fill that gap, with the goal of improving the interpretability of the model by capturing processes in the brain across groups of neurons as latent trajectories extracted from spike-train data and estimating the mutual information be- tween latent trajectories X2, XN conditioned on a behavioral variable Y, as well as the information interaction T(X1 X2 :.... XN t Y) measuring statistical dependence in X.,...... .N captured by Y.
We also aimed to extend Copula- GP with a method of ensemble model selection through the use of bagging, wherein we aggregate copula estimations of individual Copula- GP estimators fit on separate samples of data, validating the ensemble method via entropy estimation of known distribution and culminating in usage of bagged Copula-GP on aforementioned real life neuronal data, To summarize, the contributions of this paper include:
Validation of combining novel Copula-GP with GPFA (the combination of which we dub Copula- GPFA as shorthand), with the expected benefits of robustness to dimensionality and efficient, de- tailed analysis of neuronal dependence.
Justification of possible improvement to the Copula- algorithm through the addition of bagging, confirmed via validation tests on synthetic data sampled from known copula distributions.
Confirming dependencics between neuronal data from the visual cortex on pupil dilation through the use of bagged CoPAla-GPFA.
The code repository for work shown can be found at: http:://ggthub.cco//maaeee/////p
2 Neuroscience Preliminaries
In general, if system of neurons encodes the information that goes into behaviors cognition then we expect uncertainty in neuron spikes to be captured in behaviors during recorded spike-trains, or in other terms significant quantity of mutual information and/or /rr information interaction between spike signals and behaviors (Gerstner et al.) 2014). By extension, if we are able to accurately estimate mutual information then we should thusly be able to map systems of neurons to behaviors and cognition.
As stated in the introduction, computational neuroscience aims answer questions surrounding the encoding of cognition and behavior in parts of the brain. Hurdles in answering such questions in the past stemmed from the difficulty of recording groups of neurons' behavior with certainty (Kudryashova et al. 2022), however we now have the technology needed for accurate single-neuron resolution spike train trial recordings in the form of Neuropizel probes (Steinmetz el al., 2021). As such, the hurdle of neuron-behhavior-tiinne mutual information analysis is no longer centered on lack of data, moving instead to the need for accurate and efficient computational and statistical mutual information estimation techniques, as well as required robustness to dimensionality when a large number of individual neurons and behavioral variables are analysed
(Kudryashova et al.] 2022 Onken & Panzeri 2016 Vu et al. 2018
2

---

Published in Transactions on Machine Learning Research (12/2024)
3
Related Work
There are many choice methods in the realm of neuronal mutual information estimation. Popular novel ones include Bias-Improved Sraskov-Stöbbauer tion Neuronal Estimator (MINE) (Belghazi et al. 2021), both of which are effective non-parametric methods. The Copula-GP method implemented the core paper by (Kudryashova et al. 2022 is a parametric estima- tor found to generally outperform MINE and BI-KSG in mutual estimation accuracy on highly-dimensions data (2 10 dimensions), and such is our chosen method of neuronal mutual information analysis. Tt is method based copulas, multivariate distribution with uniform marginals representative of cumula- tive distributions and history of applications ranging from scientific analysis of dependency (as they are generally used in computational neuroscience (Jenison & Reale, 2004)) to the more purely analytical and predictive (i.e. predicting the outcome of an insurance claim (Hu & O'Hagan 2021 We delve more into what copulas are and how they're re useful in section 4
(BI-KSG) (Gao et al. 2018 and the Mutual Informa-
Copulas have been used in
estimation and dependency analysis in a variety of in vivo and in silicon neuronal data, including spike data, 2-photon calcium imaging, and multi-modas neuronal datasets, and have been successfully used for both mutual information estimation and dependency shape analysis in
neuronal data [Onken & Panzeri 2016 Jenison & Reale 2004
Shahbaba el. al. 2014 Berkes et al. 2008).
However, recent advances in GPU accelerated processing power (through the use of popular libraries such as PyTorch (Paszke et al. 2019 and GPyTorch (Gardner et al. 2018 have allowed for further refinement
of copula implementations and expanded the number of practical use cases, allowing for the implementation
of such method for high-dimensional neuronal data, like that of Copula-GP (Kudryashova et al., 2022)
4 Copulas: An Introduction
Copulas are multivariate distributions with uniform marginal distributions, which themselves usually rep- resent marginal variables' cumulative distributions. A d-dimensional copula function C(u,, U2, U2) [0, 1] + [] defined as a cumullatie distribution function (CDF) of a vector on the unit hyper-cube j0, 1] with uniform marginals 1I.11 (Kudryashova et al. 2022 Safaai et al. 2018). Variables of non-uniform distribution are "attached" to the marginals through the use of CDI' (which by definition map to the inter- val [0, 1), creating a copula as a joint CDF (Nelsen 2006). Sklar's 's theorem states that for a d-dimensional random vector = [X1, X2, ... ,x,} and it's CDF Fx with marginals F1, F2,... Fd, there exists a copula C' such that Vx E Rd, x ~~2 X,
Fx (11,02, 02, (2,..,,,)
m.) )
24 E.
(1)
(Sffaai et al. 2018 Van Vliet 2023) What is most valuable about this construction is the allowance for the random variables X2 Xa to take any distribution, allowing for mcaningful dependency and mutual information analysis between variables of different probabilistic distributions (Salaai et al., 2018). We can also condition the copula on some continuous variable y; allowing for the parametrization of the copula for variable like time, phase, other marginals, etc (Kudryashova et al. 2022):
Fx 82 ..... =alv)= ((F(()))))) Fa(a..,.,,, Fa(zall)v)) T. R.
(2)
To scale against higher dimensions, copulas can lake a vine copula construction, which factorize multivariate distributions into conditional distributions modelled as singular bivariate copulas (Onken & Panzeri 2016) Mitskopoulos et al. 2022). The canonical vine-construction factorises a multivariate distribution in this was
d-1d-j
( ... xa) = II (ark
II 1=1=1
1I, (a
(3)
A divide-and-conquer approach such as this allows us to mitigate the curse of dimensionality given suffi- cient performance in the bivariate case Mitskopoulos et al. (202)); Kudryashova et al.) (022). For flexible paramerization of a copula, Gaussian Process (GP) can be used for parametrization. GP essence attempts to model relationship between random variables within stochastic process such that all collections of those
3

---

Published in Transactions on Machine Learning Research (12/2024)
variables possess gaussian distribution, in effect modelling distribution of functions from which one sample Rasmussen & Williams (2005). GP can be used to estimate distribution f from which function f ~ f describing parameterization over domain [0, 1] is sampled, which can then be used in conjunction with GPLink functions to parameterize copula variants (sce table 2 in the appondix for precise definitions of GPLinks) (Kudryashova et al. 2022 Hernández- Lobato et al. 2013). To express this mathernttically, we redefine individual copulas to take the approximate form
(x ... (xal)""",,////) (zF Pa(22)02))) ... Pa(=aala(("))))0)) (1))))
(4)
where 01,... Oz, 0' T dom(y) are functions mapping some continuous variable y' onto the domain of the conditioning parameter y representing the relationships between F1y.-";2d and y', often called GPLónk functions (Kudryashova et al. 2022 Schulz et al. 2018).
Finally, we define copula entropy as
Hc(X)= -- c(u) log c(u)du,
(5)
where c(u) represents the probability density of the copula at u
JCC
c(u)
(6)
Du, Guz .elg....uu. Ous
(Ma & Sun 2011 Jenison & Reale 2004) As proven by Ma & Sun 2011). copula entropy equivalent the negative mutual information of the marginals, i.c. T(X) -H (X). In other words, through accurate
estimation of larget distribution's copula we are able to accurately estimate the mutual information between that distributions' marginal variables.
5 Copula-GP
The Copula-GP package was implemented in Python and deployed in Kudryashova et al. (2022). and is one of the few (if not only) practical implementations of GP-treated copula vine model designed with GPU based acceleration of computation via Pytorch (Paskke al. 2019). is also the first of such models designed with use on neuronal data in mind, and as such literature the uses of GP vine copulas modelling neuronal data is heavily limited, with the main source for such claims being the core paper. In simplest terms, the package serves as framework for GP copula model parameter estimation, copula model selection, and model deployment as a vine, with additional built in visualization and entropy extraction utility. The core paper also utilizes comparison to other popular methods of neuronal mutual information extraction, yielding superior results in highly-dimensional data. Any subsequent claim surrounding the Copula- GP package without citation attached in this section comes courtesy of (Kudryashova al.) 2022).
The vine construction implemented in Copula-GP used copula building blocks from five distinct families: independence, Gumbel, Gaussian, Frank, and Clayton copulas, with the independence copula preferred for independent variables (see figure G The factorization of the vine is that of the C-vine described by equation g with accommodations made for putaaeramttriiniinnn. To fully capture the tail dependencies and negative correlations in relationships between marginal distributions, mixed copulas were utilized. In the core paper, these defined as
K Cmiied ( = Σ oy (YC,((11))))),
(7)
where K is the number of elements, oj is the concentration of the jth copula (cz), and 0 is the jth copula's 's parameter GPLink. The GPLink for each copula determined by its copula variant, as shown in figure 2
with 0 being defined by GPLink(J), where f is sampled from 0 ~ Niu, K2 (X, X)) (the choice of GPLink
4

---

Published in Transactions on Machine Learning Research (12/2024)
depends on the kind of copula; see 2b. GP is also utilized to parameterize the concentrations Ø which are defined as
+j - (1-4)
(4)
II mm Ln o [m +
+ Φ- M - m ( LM 0,
(8)
where & is the CDF of a standard normal distribution and fm is sampled from fm ~ (i K (Y,)). This gives us 2M - sets hyper parameters to estimate, {A}M kernel hyperparameters for each GPLink 0 and (A/as- kernel hyperparamcters for cach concentration function Φ, estimated via the methods described in section E.
6 Entropy and Mutual Information Estimation
As Copula- GPpula-GPaal GP models copula distribution, possible states of variables can be sampled from it. As such, the joint mutual information r(x1 X² :...: ..........
= :)) distribution X = {X1, _ can be derived from C(X(Y) (the C-vine copula Copula-GP estimates). Copula-GP possesses built-in capabilities for computing mutual information directly via the above equation, however doing so requires computation of nested integrals and such is computationally intensive; in testing, when attempting to extract the mutual information of complex distributions we found our machine quickly exhausted GPU memory when utilizing single NVidia 2080Ti. We choose to find the MC estimate (Robert & Casella 1999) of Hx(X) (which is also utilized in the core paper), which invol ves sampling from the distribution C(X (0:) with parameterization in S random values ÿJ ~ (((0 1) and estimating the conditional entropy as the mcan
(x)
Hc(X(y) between parameter Y and multivariate
Ho(X = f( Ho(Xi R
5
Hc(X
(9)
dom(y)
Similarly, we estimate the unconditional entropy as the mean entropy of the system found when fed true values of the paramcterizing variable y: 10, i.c
1 S Σ Hc(()yn).
Hc(X) = (Hy(Hc(XlV) ~
(10)
Due to the law of large numbers, as N goes to infinity both approximations become more accurate. To find the point-wise copula entropy values ()(()),), we utilize Copula-GiP GP's implemented entropy extraction methods
7 Extension of Copula-GP: GP: Bagging
Bagging is bootstrapping technique wherein multiple "weak" models are trained and their outcomes ag- gregated in some way to create a "stronger" model, typically reducing variance and increasing model bias. We propose a formal justification of copula bagging, taking inspiration from the "Random Forest" algorithm for bagged regression (Breiman 20010; we assume the existence of some true copula function C describing a multivariate distribution X. We take N samples ( , and get a copula ((n) representing the best copula fit possible on an individual sample Sy, We then come final estimate via mean aggregation of our estimates. Suppose ((n) C-vine estimate, and let city) represent the i-th copula in the j-th layer of the n-th estimate. We can create a mcan mixture copula Ci representing the final estimate for the corresponding true copula Cr,j as Č₁, = N Σ (n) As N + oo, by the law of large numbers (and the fact that copulas are distributions) we find Ct + Cijz and so Ô - C. In other words, by fitting copulas on individual spike-train trials and taking the aggregate, we estimate a copula that more closely aligns with the true probabilistic relationship than those aggregated.
7.1 Implemented Weighted Aggregation Methods
A simple unweighted mean is appropriate in the infinite case. In practice however, we often find that a naive average insufficient when examining finite number of samples; suppose we have finite N uniform length
5

---

Published in Transactions on Machine Learning Research (12/2024)
samples Sn, and thus have set of observations S = C.. If the n-th sample is an outlier sample and so Cn does not resemble the true copula Cn, we naturally should account for that. As such, in addition to a naive average, we may choose to weigh models based on various Bayesian criterion. We utilized a Bayesian model aggregation approach utilized by S. Hu et. al (Hu
fitting corresponding copula estimates
& O'Hagan 2021), which aggregated copulas via Bayesian Information Crilerion (BIC). The BIC value for the n-th estimated copula given as
BIC, = 2log ((C,|X,)) + Pn log X
(11)
where log L(C,,) X; @) is the log-likelihood of copula Č₂ under observations X and estimated model parame- ters , and Ps is the total number of estimated parameters of Čn (in the case of Copula-GP, each bivariate copula. estimated is a mixture copula and so possesses mixing and dependence parameters for each copula mixed; if the mixture is singleton mixture, we only count dependence parameter. Note independence cop- ulas have no dependence paramcter). BIC essentially rewards model log-likelihood with penalty in number of parameters scaling with sample-size (Schwarz 1978 Hu & O'Hagan 2021), with small BIC values are being preferred (as negative as possible). We arrive at weights for the n-th model as
exp( exp(-
W = Σm 2xp((2)
(12)
(Hu & O'Hagan 2021) We also utilize the Akaike Information Criterion (AIC) for comparison, which is given for the n-th copula estimation as AIC, = -2 log [(X]Cn, Q) + 2pk and extract weights AI, accordingly via equation 12 replacing BIC, with AIC,. The AIC is similar to the BIC in that it rewards log-likelihood but has a more relaxed penalty in number of parameters. If a more complex copula matches the true distribution, this can lead to possible improvements in aggregation accuracy, but can also lead to over-fitting via wwerpparameterization (Zhou 2021 (as we are already utilising mixtures of copulas with a high ceiling in number of parameters, we expect the latter to be true). We additionally include dynamic weighting of estimates via calculating the above information criteria and their respective weights point-wise as opposed to setting weights to be constant over the observations X; we call point-wise aggregation dynamic bagging and constant weighting static bagging. We validate our bagging methods with different aggregation methods on data samples from randomly generated bivariate copulas in section g
To aggregate C-vines, we bag one layer at a time; we first bag copulas in the current layer via our selected mothod and gather cumulative conditional probabilitics as pscudo-observations via the bagged copulas' conditional cumulative distribution functions (cedfs). We then utilize them for BIC, AIC, and R² calculations in the next layers' bagging process, and doing so effectively propagate the previous layers' weightings to the next layer. As subsequent bagged copulas will have worse BIC, AIC, and R caleulations if their corresponding C-vines' previous layers were found to be bad fits, we in effect prioritize relationships with weaker conditioning using this method.
Copula-GP estimates a parametric C-vine along given 1-D parameterization input X, modeling each individual copula as a mixture and getting estimated dependence parameterizations (X), and estimated mixing parameterizations (Y), for each of i mixture copulas (with each mixture copula being a mix of j(o) copula variants). As such, our actual implemented bagging procedure to aggregate along these estimated parameters for each unique copula variant in each mixture (we count each rolation variant of a algorithm Clayton or Gumbel copula as unique; this algorithm is outlined more formally in the appendix,
D. In addition, implemented model selection for each C-vine we aggregate, allowing for the bagged vine to contain copula mixtures possibly unexamined in the heuristic based model selection process. The C-vine aggregation process requires aggregation of copulas by layer and storing the cumulative distributions of each copula via the bagged copula ccdfs for weight extraction in the next layer.
8 Canssian-Process Factor Analysis
A very common tochnique in the realm of dimension reduction for highly dimensional neuronal data is Cilussian-Proces Factor Analysis (GPFA), which factors a matrix of observations X into a product of a loading matrix V and scoring matrix θ as X = VO (the Factor Analysis step) producing a low-dimensional
6

---

Published in Transactions on Machine Learning Research (12/2024)
group of trajectories, and then smooths said data by means of fitting a GP (the Gaussion Process step) Pirs
& Strumbelj 2022). This operates under the assumption of a linear relationship between observations X and latent trajectories θ (utilizing V/ as a transformation matrix) (Yu et al. 2009). and models the GP with bias as
x= V ((,
(13)
with parameters (d, V The Copula-GPFA process consists of GPFA applied to neuronal data to extract latent Irajectory estimations X1, X2,..., followed by fitting a Copula-GP C-vine estimator on neuronal trajectories. Doing so, we are able to efficiently extract mutual information estimates describing not just the information interaction between the examined part of the brain and some parametcrizing variable, but the mutual information between the brain processes driving neurons instead of individual neurons themselves (Yu et al. 2009: Ma & Sun, 2011 Kudryashova et al. 2022), enhancing the interpretability of individual copulas estimated. In addition, by utilizing GPFA to reduce a large number of neurons to a more manageable amount, we reduce the computational intensity of mutual information analysis; we extract 13 trajectories from 155 neurons then from the naive estimate for the complexity training Copula-Gi estimators (see core paper [Kudryashova et al. 2022 we can come an upper-boune bound of 150² 13a 133 times faster Copula-GP training.
9 Bagged Copula-GP Validation Tests
To validate our bagging methodology, we tested various Copula-GP aggregation methods alongside baseline unbagged Copula-GP on generated bivariate copulas (the "true" copula). C-vines made up of bivariate copula building blocks, with C-vine model selection essentially being made up of consecutive bivariate copula selections. As such, validation alongside bascline Copula-GPP on bivariate copulas alone is indicative of improvements in C-vine entropy estimation via bagging methods described.
For all validation tests in this section, the task was to accurately replicate the true copula's parametric entropy and dependency shape, parameterized on some variable æ. For validation, the weighted aggregation methods utilized were 1. Copulas weighted point-wise dynamically by BIC / ATC on input, 2. Copulas weighted statically by BIC AIC input, and 3. naive average of copulas. In addition, we also examine copula entropy point-wise root mean squared error (RMSE) of baseline and bagging methods.
To generate the data, we drew 10000 samples from the true copula, splitting samples and their corresponding parameterizations into train and test sets (80:20 split). Validation tests and 2 utilized parametcrization on a normalll distributed variable x ~ (0.55, 0.2,, restricted to the interval [0, 1]. Validation test 3 utilized a parameterization in a time-based variable scaling linearly from to Copula training and model selection utilized train set samples, with accuracy validation utilizing metrics described on test set samples. All validations utilized the heuristic algorithm for individual estimator model selection, and for bagged imaatioo 4 estimators underwent individual training and model selection routines. We include random seed, module version, and hardware information for test reproducability purposes in the appendix (see G
Table 1: True and predicted mean copula entropies He (with 95% confidence interval included in parentheses) and point-wisc RMSE of validation tests. Closest to actual best scores bold. We find that BIC Dynamically bagged Copula-Gr GP consistently outperforms alternative methods in accuracy.
Model
Test He
RMSE
Test H
RMSE
Test 3 H
RMSE
True Copula
-0.1440 (0.0154)
-0.6864 (0.1148)
-0.6864 (0.1148)
Baseline
-0.0000 (0.0000)
0.1442
-0.6680 (0.0688)
0.0694
-0.4556 (0.2038)
0.2587
BIC Dynamic
-0.0802 (0.0147)
0.0647
-0.6764 (0.0847)
0.0719
-0.5628
(0.0754)
0.1388
BIC Static
-0.0791 (0.0117)
0.0656
-0.6725 (0.0793)
0.0713
-0.5589
(0.0650)
0.1418
AIC Dynamic
-0.0796 (0.0128)
0.0652
-0.6664 (0.0870)
0.0748
-0.5382
(0.0813)
0.1615
AIC Static
-0.0790 (0.0123)
0.0658
-0.6662
(0.0921)
0.0658
-0.5389
(0.0703)
0.1596
naive average
-0.0798 (0.0133)
0.0651
-0.6738 (0.0778)
0.0651
-0.5580
(0.0722)
0.1420

---

Published in Transactions on Machine Learning Research (12/2024)
Low-Entropy Complex Mix Validation Test
Baseline Estimation
Dyn. BIC Weighted Bagged Estimation
True Copula
1.0
1.0
1.0
0.8
0.8
0.3
0.6
0.6
0.4
0.4
0.4
0.2
0.2
0.2
0.0
0.0
0.0
0.00
0.25
0.50
0.75
1.00
0.00
0.25
0.50
0.75
1.00
0.00
0.25
0.50
0.75
1.00
Independence
Indep. ClaytonO +90
Indep. Clayton 90"
+ +180*+270" Gaussian
Gaussian Gumbel 0
Figure 1: Average predicted copulas of BIC dynamically bagged and baseline Copula-GP on low-entropy copula validation test set, with true copula for comparison. Contour lines of copulas are visible. Note how the distribution "tricks" the bascline estimator into picking the independence copula, but the bagged estimate is able capture major tail dependencics in the true distribution (circled in red).
9.1 Results
For the first validation, we expected our bagging methods to perform well. Baseline Copula-GP when facing such a copula can often "give up" early into heuristic model selection and select independence if it finds non-independence copulas possess too low WATC. We expect that bagging able to make up for this flaw in baseline heuristic model selection via it's bias-variance trade-off, and thus pick up tail dependencies in the model. In table we observe that bagging methods all performed better than baseline (which selected independence), with BIC dynamic bagging possessing the best model accuracy and RMSE. Finally, in figure we find that the BIC dynamically bagged estimated copula captures major tail distributions in the true copula, whereas the baseline estimation of independence does not.
High-Entropy Simple Validation Test Dyn. Weighted Bagged Estimation
Baseline Estimation
True Copula
0.8
0.8
0.8
0.6
0.6
0.6
0.4
0.4
0.4
0.2
0.2
0.2
0.25
0.75
0.25
0.75
Clayton90 Gumbel90"
Gaussian Gumbel90* +270" Clayton90*
Gumbel90' Clayton90'
Figure 2: Average predicted copulas of BIC dynamically bagged and baseline Copula-GP high-entropy copula validation test set, with true copula for comparison. contour lines of copulas are visible. We find that bagging hods can result in incorrect selection of copula variant leading to possible misrepresentation of dependency shape, however the bagged estimation outperforms baseline mean copula entropy estimated.
For the second validation test, we expected baseline perform well and were more interested in if the bias induced by bagging result in worse performance distributions the baseline is well suited for. The true copula here is high-entropy copula with simple shape, something that in the core paper (Kudryashova ct al. 2022) was found to be well suited to bascline Copula-GP as a method. The BIC dynamically weighted estimate had the best. estimated mean copula entropy, with the naive average estimate possessing superior RMSE. In figure we find that the BIC dynamically weighted estimation included a gaussian copula and as such makes incorrect assumption surrounding the dependency shape of the true copula. Despite this, the majority of bagged copula estimates still result in similar if not slightly better performance in copula entropy prediction.

---

Published in Transactions on Machine Learning Research (12/2024)
High-Entropy Transformed Copula Validation Test
Baseline Estimation
Dyn. BIC Weighted Bagged Estimation
True Copula
1.0
1.0
1.0
0.8
0.8
0.8
0.6
0.6
0.6
0.4
0.4
0.4
0.2
0.2
0.2
0.0
0.0
0.00
0.25
0.50
0.75
1.00
0.00
0.25
0.50 0.75 1.00
0.00 0,25 0.50 0.75 1.00
Gaussian Clayton90* +2700
Gaussian Gumbel90s+270"
Gumbel90* +2700" Clayton90*
Clayton90 +2700
Figure 3: Average predicted copulas of BIC dynamically bagged and baseline Copula-GP on transformed high-entropy copula validation test set, with true copula for comparison. contour lines copulas are visible. We find that the baseline estimation possesses less of resemblance to the true copula than the baggod estimation.
The final validation was on a copula parameterized time. In the core paper, it was found that "transformed" cases with dependency shape changing function in time can be more challenging for Copula-GP to predict the entropy of accurately. As such, utilize the copula in validation test and change its parameterization to be in time t scaling linearly from O to 1. We maintain continuity when splitting into train and test set, and as such the train set consists of parameterization parameterization
f = 0.0 to 0.8, and the test set consists of t = 0.8 to 1.0. We hope that bagging distribution estimations made by cstimators trained on continuous subsets of training data can induce some robust ness time-based transformation true dependency shape. Our hopes were confirmed, we find in table ! that bagging methods all find better accuracy than bascline. BIC dynamically weighted bagging results in the best mcan copula entropy estimate and RMSE. From figure 3 we see that the baseline estimated copula possesses much less of resemblance to the true copula than the bagged estimate.
9.2 Validation Discussion
In general, we found that bagged copula estimations utilizing 4 estimators possesses similar or better accuracy in copula entropy estimation than baseline. BIC dynamically aggregated copulas possessed the best accuracy for all validation tests conducted, often with the best point-wise R.M... We also found that naive average of copulas estimated performed surprisingly well, often times out performing other bagging methods. From validation test 1, we see that bagged estimates can catch dependencies in distributions that the baseline perceives as independent. In C-vine model selection, as copulas gain more conditioning marginal variables' relationships tend to more closely resemble conditional independence (Kudryashova et al. (022). As such, we may conclude that bagged estimator might be more likely pick up tail dependencies in conditioned relationships where the baseline estimator might assume independence. From validation test 3 we find that the bagged estimator is more robust to transformations in dependency shape over time, and as such may yield better accuracy when predicting marginal variables' dependency shape for time-scrics. Finally, we note that the BIC dynamically bagged estimator does not over-estimate copula entropy in validation.
10 Exploration of Experimental Neuronal Dataset
For this section, the objective was to utilize Copula-GPFA to confirm high statistical dependence between the visual cortex and pupil dilation. The main data set we wished to explore is the Visual Coding: Newropizels dataset, a publicly available dataset containing spike-signals recorded from the visual cortex of live mice at single neuron spatial resolution utilizing novel Neuropixel probes, a high-fidelity and -resolution brain probe developed 2021 (Stcinmetz al. 2021). Due to the volume of in vivo data recorded, as well the quality metrics included in the data, it an invaluable tool for analysis of underlying processes within the visual cortex. All technical claims surrounding the data in question is sourced from the technical white
9

---

Published in Transactions on Machine Learning Research (12/2024)
paper (sce Dj. For more information surrounding data production and specificities not covered in the body of this paper, i. specifics on how probes were used and how visual stimuli were produced, one should see the technical white paper. The spike data was collected via insertion of Newropizel probes into the visual cortex of live mice. Probes record spike signals from points in space at single ncuron resolution. These points are called "units," " and are where neurons might be. Units have quality metrics attached, which can be used to filter for units of better quality (i.e. units with less noise, units which more certainly record spikes from only one neuron, etc). We chose units via the cleaning processes described in D utilizing these quality metrics, arriving at dimensionality of 178.
Pupil dilation has been found to have close ties to processes within the visual cortex (Franke et al., 2022 Bombeke et al.) 2016 Larsen et al.) (018)) As this data is directly from the visual cortex, we expected a significant statistical dependence between pupil dilation and visual cortex trajectories. As such, the motivation behind utilizing this dataset is as an experimental method of model validation. We hypothesisod significant information interaction extracted via Copula-GPFA, thereby effectively confirming that the visual cortex of the brain is linked to pupil dilation and agreeing with prior results in the literature (Bombeke et al.
2016
Ganea et al. 2020).
10.1
Runtime Breakdown
Runtime of Copula-GPFA can be broken down into three steps: transformation into latent trajectories via GPFA, fitting of multivariate distribution via Copula-GP on either part of whole data, and aggregation via BIC dynamic bagging. GPFA was fit on 100 consecutive spike trains during Drifting Gratings stimulus presentations, each of which have uniform temporal length of 2s followed by a is interval during which no stimulus is presented. We fit the model for the entirety of the dataset with the dimensionality of the elbow of the log-likelihood plot (figure 9). reducing the original dimensionality of 178 down to a latent dimensionality of 13. After GPFA application, little interim processing is needed: we removed drift and concatenated trials as per the process described in section F.2 to produce remi-continuous data and mapped onto the domain (0, 1) via additional interim steps described in scction F.2
We utilized Copula-GP to estimate C-vines over the concatenation of GPFA-treated trials, possessing 12 layers (78 copulas total) with parameterization in normalized pupil dilation (as per steps described in D testing both bascline unbagged and bagged Copula-GP. For bascline, Copula-GP was fit over the full semi- continuous data of 80 trials. For bagged Copula-GP, 4 individual models were fit over semi-continuous subsets of 20 trials each.
For model ensembling, we utilized the BIC meeamicaall-weighhae aggregation method. We estimated para- metric negalive copula entropy (XIY) over time utilizing baseline and bagged estimates for 15 continu- ous trials, expecting either similar or lesser (more negative) copula entropy estimates of bagged Copula-GP as baseline, as was observed in validation.
10.2 Results
We can see the outcome of - (X) and interaction information T(X +- Y) calculations alongside normalized pupil dilation in figure 4 Bagged Copula- GP estimated a copula with mean entropy (X) = 13.4027
bits (95% CI of 0.7069) and mean - (X)(X = 7.3323 bits (95% CI of 5.9693). Copula entropy appears to be highly correlated with pupil dilation, with decreases in copula entropy occurring with decreases pupil dilation. In addition, the information interaction stays firmly negatine, implying the multivariate distribution of neuronal trajectories possesses statistical dependence on pupil dilation (thereby confirming our hypothesis).
From figure 5 we can see that like the bagged estimate the bascline negative copula entropy appears to be dependent on pupil dilation, with large dips occurring with pupil dilation. We found a mean negative copula entropy (X) estimation of 13.4745 bits (95% CI of 3.0720) and a mean negative conditional copula entropy (X)) empirical estimation of 7.2221 bits (95% CI of 6.1190), with X being the neuronal trajectories and Y being pupil dilation. However, we find the bagged estimate is less sensitive to small changes in pupil dilation than the baseline estimate. In addition, the bagged entropy estimate appears
10

---

Published in Transactions on Machine Learning Research (12/2024)
to steadily hover around the mcan negative entropy found for most observed time-buckets where the unbagged estimate appears almost linear in pupil dilation (see figure These differences aside, the bagged and baseline mean copula entropy and information interaction estimates were close, implying that the bagged estimator does not inherently out- or under-perform baseline unbagged Copula-GPFA in the case of the in vivo data.
Bagged Est. -Hc(X) (bits) w/ Pupil
Bagged Est. 1[x + Y) (bits) w/ Pupil Dilation over 15 Continuous Trials
Dilation over 15 Continuous Trials
14.0
-4.0
IEX - Y0
0.625
Pupil Dilation (Norm.)
IX Y7
0.600
13.0
0.575
12.5
M
0.550
-5.5
12.0
0.525
-6.0
11.5
(Me(0)
0.500
-Mc00)
-6.5
0.475
200
400
600
800
1000
1200
1400
200
400
600
800
1000
1200
1400
time bucket
time bucket
Figure 4: Pupil dilation (in red), copula entropy of neuronal trajectories parametric in pupil dilation utilizing bagged estimation (in blue), and corresponding information interaction (in green) through time utilizing estimated copula entropy. Copula entropy was estimated via estimating a distribution with BIC dynamically bagged Copula-GP.
Bagged vs. Unbagged Est. -Xc(X) (bits)
Bagged vs. Unbagged Est. 1(X < Y) (bits)
15
14
1(X Y
12:
-5
Bagged (Xc(X)
Bagged -HcO0
Unbagged -H²(X)
10
Unbagged -Hc(X)
200
400
600
800
1000
1200
1400
200
400
600
800
1000
1200
1400
time bucket
time bucket
Figure 5: BIC dynamically bagged vs unbagged Copula- GP estimated entropies and information interactions. The BIC dynamically bagged entropies curve seems to hover steadily around mcan, and does not appear to be as volatile as unbagged entropies.
11 Discussion
Copula-GPFA able to efficiently provide accurate and meaningful dependency and infor- mation quantification analysis and results. We validated Copula-GPFA as mutual information and information interaction extraction method on in vivo data extracted from mice, confirming significantly nega- tive information interaction between neuronal trajectorics and pupil dilation and agrecing with prior findings linking the visual cortex with pupil dilation (Bombeke et al.) 2016 Ganea et al. 2020 Zylberberg el al. findings, allowing for visualization of both neuronal trajectories and information quantities over time. We also justify a significant reduction computational complexity via GPFA utilized with Copula-GP, using complexity estimates made in the core paper.
2012). In addition, Copula-GPP GPFA also allows for easy interpretability
11

---

Published in Transactions on Machine Learning Research (12/2024)
Dynamically bagging Copula-GP estimates over subsets of data provides superior estimations over baseline in tail-dependency identification and estimated copula entropy accuracy in the bivariate case. We also introduced copula bagging as a means to possibly improve copula estimation, complete with both formalistic and practical justifications, and implemented bagging as an extension to the Copula-GP code-base (locally). In addition, we utilized the bagged exlension on the experimental in vioo data, finding results in agreement with baseline estimation of mean copula entropy with more stability over time, however whether this stability is indicative of the true distribution or not goes unconfirmed. As Copula-GPFA is not inherently method restricted to the area of neuroscience, it is entirely plausible Copula-GPFA is able to be meaningfully applied to other fields where dependency analysis is useful (i.e. bio-informatics or quantitative finance).
While we did not test robustness of the bagged extension to higher dimensionality, as C-vine selection within the Copula-GP algorithm consists of individual selections of bivariate building blocks we hypothesize that benefits found in the bivariate case propagate into cases of higher dimensions. Our implementation of bagging also follows in the footsteps of the formalistic approach of the original Copula-GP implementation, and such maintains flexibility in use cases outside the realm of spike train data and neuroscience. We also acknowledge that the ability to bag copula estimations allows for the aggregation of Copula-GP estimators to be fit on individual trials of data, thereby resolving possible violations of local smoothness and continuity assumptions caused by concatenation of GPFA treated spike train trials (due to time-to-fit. concerns and limitations in compute-power, this goes unutilized as of now). Finally, we note our bagging implementation can be used without heuristic model selection, and one could utilize our implementation to bag singleton mixtures (mixture copulas of only single variant) as an alternative efficient. means of model selection. Finally, we note that while we did not thoroughly investigate individual bivariate copulas within the C-vines estimated, each bivariale copula ilself represents dependency relationship. As such, there is room for more analysis.
11.1 Further Validation, Optimization, and Use of Bagging Methods
Our validations find that weighted copula estimate aggregation can yield effective improvements in copula tail dependency identification and copula entropy accuracy. The obvious extension to validations made would of course be to confirm if such improvements firmly carry over to higher dimensions, however because Copula-GP implementation and the structure of the C-vine we believe that benefits found the bivariate case will propagate into higher dimensions. In addition, further refinement of aggregation methods may be considered, i.c. clustering methods such as k-means Use of parallel processing in training bagged estimators concurrentyy might also be a natural optimization for bagged Copula-GP; one could even bag bivariate copulas immediately after they've completed model selection on their subsets of data during vine training.
way to select copula estimates for datasets.
That being said, benefits found in the bivariate case are still entirely applicable. Bivariate copulas are still effectively used in fields from computational finance (Cherubini et al.] [001]] to bio-informatics (Ray et al., 2020). As such, the improvements made to Copula-GP in bivariate copula estimation via the addition of dynamically weighted aggregation methods may be used to robustify bivariate copula selection effectiveness in such use cases.
Finally, we believe there is a gap within the literature for a more comprehensive review of mutual information estimation. While Copula-GP was tested against heh commonly used methods in the core paper (KSG by Kraskov et al. (2004), BI-KSG by Gao et al. (2016), and MINE by Belghazi et al. (2021)) 1), other methods have recently been produced with robustness to high- digh-iimensionlltty. In particular, several recent novel Density Ratio Estimation (DRF) based methods have been produced, such as TRE by Rhodes et al. (222)) DRF-
by Choi et al. (2022) and MDRE by Srivastava et al. (2023)
References
Mário Antunes, Diogo Gomes, and Rui L. Aguiar. Knee/ Elbow Estimation Based on First Derivative
Threshold. In 2018 IEFF Fourth International Conference on Big Data Computing Service and Ap-
12

---

Published in Transactions on Machine Learning Research (12/2024)
plications (BigDataServicee) pp. 237-240, March 2018. doi: 10.110//BsDDaaaSrrrce..0000000 URL
https://iseexpporrs ieee. org/ http:///pppppppeerrrrrr////////p /ttp:///aaccc//cccctt//////7/7
Bence Bagi, Michael Brecht, and Juan Ignacio Sanguinetti-Shheck Unsupervised discovery of behaviorally
relevant brain states rats playing hido-and-sock. Current Biology, 32(12):2640-26533 June 2022.
htt.://wwwweeeeeeee...omoooeeeee ISSN 09609822. doi: 10.1016/j.cub.20222 (tt.ccub.222..4.0000 URL http com/retrieve
pi1/80960962220070777
Mohamed Ishmael Belghazi, Aristide Baratin, Sai Rajeswar, Sherjil Ozair, Yoshua Bengio, Aaron Courville, and R. Devon Hjelm. MINE: Mutual Information Neural Estimation, August 2021. URL http://arxiv.
org/ hrg/abs/18011.. 04062 arXiv:1801.0406.
P. Berkes, F. Wood,
and J. Pillow.
Characterizing neural dependencies
with copula models.
December 2008
URL https://wwrr www.sntcscccchhhlrrr........
Gérard Biau and Erwan Scornet. A random forest guided tour. TEST, 25(2):197 227, June 2016. ISSN 1863-8260. doi: 10.1077/111488-01--0453 URL
Inga Blundell, Romain Brette, Thomas A. Cleland, Thomas G. Close, Daniel Coca, Andrew P. Davison, Sandra Diaz-Pier, Carlos Fernandez Musoles, Padraig Gleeson, Dan F. M. Goodman, Michael Hines, Michael W. Hopkins, Pramod Kumbhar, David R. Lester, Bóris Marin, Abigail Morrison, Eric Müller, Thomas Nowotny, Alexander Peyser, Dimitri Plotnikov, Paul Richmond, Andrew Rowley, Bernhard Rumpe, Marcel Stimberg, Alan Stokes, Adam Tomkins, Guido Trensch, Marmaduke Woodman, and Jochen Martin Eppler. Code Generation in Computational Neuroscience: Review of Tools and Tech- niques. Prontiers in Neuroinformatics, 12, November 2018. ISSN 1662-5196. doi: (0.33899/minf.201... 00068. URL attps://wxit www.:/footieisiinocmm http://w...................ooota
full Publisher: Frontiers.
Klaas Bombeke, Wout Duthoo, Sven C. Mueller, Jens-Max Hopf, and C. Nico Boehler. Pupil size di- rectly modulates the feedforward response in human primary visual cortex independently of attention. Newrolmage, 127:67-73, February 2016. ISSN 1053-8119. doi:
URL
1023/A:10109334443244 URL httpp//ccccccccc/cocccccccccccca
Leo Breiman. Random Forests. Machine Learning, 45() 32, October 2001. ISSN 1573-0565. doi: 10.
Umberto Cherubini, Sabrina Mulinacci, Fabio Gobbi, and Silvia Romagnoli. Dynamic Copula Methods in Finance. John Wiley & Sons, November 2011. ISBN 978-0-470-68307-1. Google-Books-ID: wZaAD-
WAAQBAJ.
Kristy Choi, Chenlin Meng, Yang Song, and Stefano Ermon. Density Ratio Estimation via Infinitesimal Classification, March 2022. URL http://arxiv.. org/abs/2111. 11010, arXiv:2111.11010.
Gabricla Czanner, Sridevi V. Sarma, Demba Ba, Uri T. Eden, Wei Wu, Emad Eskandar, Hubert H. Lim, Simona Temereanca, Wendy A. Suzuki, and Emery N. Brown. Measuring the signal-lo-noise ratio of neuron. Proceedings of the National Academy of Sciences, 112(23): June 2015. doi: 10. 1073/pnae
http:/aaa.11555551111 1505545112. URL httpp://www. pnas. /ttt:///////b///..1///////////// 1505544112 Proceedings of the National Academy of Sciences.
Publisher:
Michael Denker and Moritz Kern. Flephant 1.0.0, November 2023. URL https://zenodo: http:://zenddd.r///eerras///
org/reoords/
10107318
Daniel A. Dombock, Anton N. Khabbaz, Forrest Collman, Thomas L. Adelman, and David W. Tank. Imaging largo-scale neural activity with cellular resolution in awake, mobile mice. Newron, 56(1):43-55 -3-)7,, October 2007. ISSN 0896-6273. doi:
13

---

Published in Transactions on Machine Learning Research (12/2024)
Katrin Franke, Konstantin F. Willeke, Kayla Ponder, Mario Galdamez, Na Zhou, Taliah Muhammad, Saumil Patel, Emmanouil Froudarakis, Jacob Reimer, Fabian H. Sinz, and Andreas S. Tolias. State- dependent pupil dilation rapidly shifts visual feature selectivity. Nature, 610(7930):128 134, October 2022. ISSN 1476-4687. doi: 10. (188//415860022-052777
URL https://www nature www.coo/cca//o articcesteess/ +41586-022-00270-3 Publisher: Nature Publishing Group.
Dan Alin Ganea, Alexander Bexter, Mathias Günther, Pierre-Marie Gardères, Björn M. Kampa, and Florent Haiss. Pupillary Dilations of Mice Performing a Vibrotactile Discrimination Task Reflect Task Engagement and Response Confidence. Frontiers in Behavioral Neuroscience, 14, September 2020. ISSN 1662-5153. doi: wttp://nnehh220000000
URL http::/www..rttttrrrr.rrrrrraaaa http://ricceee/11....388//////// 2020.
00159
Publisher: Frontiers.
Weihao Gao, Sewoong Oh, and Pramod Viswanath. Demystifying Fixed k-Nearest Neighbor Information Estimators, August 2016. URL attp://arxiv.org httpaawww.rggga////...00000 wttpaaaba......000 03006 #xXiv:160403306
Weihao Gao, Sewoong Oh, and Pramod Viswanath. Demystifying Fixed k -Nearest. Neighbor Information Es-
timators. JEEE Trunsactions Information Theory, 64(8):5629 5661, August 2018. ISSN 1557-9654. doi:
10.1109/TIT.200 2807481 URL https://ioeexplorrs leee. wwgg/a/ssraatt/ocmmm/////22222 Con-
ference Name: TEEE Transactions on Information Theory.
Jacob Gardner, GeofT Pleiss, Kilian Weinberger, David Bindel, and Andrew G Wilson. GPyTorch: Black- box Matrix-Matrix Gaussian Process Inference with GPU Acceleration. In Advances in Neural Information Processing Systems, volume 31. Curran Associates, Inc., 2018. URL https: stttop://proceeedngg. wewpedddiinnnnnuiis..... nwuripurrrcccc Abstract html
Wulfram Gerstner, Werner M. Kistler, Richard Naud, and Liam Paninski. Neuronal Dynamies: From Single
Newrons to Networks and Models of Cognition. Cambridge University Press, July 2014. ISBN 978-1-107-
06083-8. Google-Books-ID, D4j2AWAAQBAJ.
Joshua T. Glaser, Ari S. Benjamin, Roozbeh Farhoodi, and Konrad P. Kording. The roles of supervised machine learning in systems neuroscience. Progress in Neurobiology, 175:126 137, April 2019. ISSN 0301- 0082. doi: 10.1016/j. wtwp://jpppeuroboo... 201...0....11.....000 URL http:://ww..ccceeeddreccc http://d........ooocccc//////aaa pi1//03000002213308666
José Miguel Hernández-TTobtto James R Lloyd, and Daniel Hernández-Lobaue Gaussian Process Condit ional Copulas with Applications to Financial Time Series. July 2013.
Sen Hu and Adrian O'Hagan. Copula Averaging for Tail Dependence in Insurance Claims Data, March 2021. URL attp:://axxle wttpi//.orrgaa////.s//// 10912 arXiv:2103. 10912 [stat].
Aisha Jangid, Laxmi Chaudhary, and Komal Sharma. Computational Neuroscience and Its Applications: A Review. In Mohammad Shorif Uddin, Avdhesh Sharma, Kusum Lata Agarwal, and Mukesh Saraswat (eds.), Intelligent Energy Management Technologies, pp. 159-169, Singapore, 2021. Springer. ISBN 158824. doi: 10.1007 070007998--811--88888
16.
Rick L. Jenison and Richard A. Reale. The Shape of Neural Dependence. Neural Computation, 16(4):665- 672, April 2004. ISSN 0899-7667. doi: 10.1162/08997668 /099099976666049226666655 URL https:///oi.orggee 1162/ 089976604322866559
James J. Jun, Nicholas A. Steinmetz, Joshua H. Siegle, Daniel J. Denman, Marius Bauza, Brian Bar- barits, Albert K. Lee, Costas A. Anastassiou, Alexandru Andrei, Çağatay Aydı,, Mladen Barbic, Tim- othy J. Blanche, Vincent Bonin, João Couto, Barundeb Dutta, Sergey L. Gratiy, Diego A. Gutnisky, Michael Häusser, Bill Karsh, Peter Ledochowitsch, Carolina Mora Lopez, Catalin Mitclut, Silke Musa, Michael Okun, Marius Pachitariu, Jan Pulzeys, P. Dylan Rich, Cyrille Rossant, Wei-lung Sun, Karel Svoboda, Matteo Carandini, Kenneth D. Harris, Christof Koch, John O'Keefe, and Timothy D. Har- ris. Fully integrated silicon probes for high-density recording of neural activity. Nature, 551(7679): 232-236 236, November 2017. ISSN 0028-0836, 1476-4687. doi: 10.1038 1.11088/natrr2466. URL https:
7/www. wwwwwentur...mmm com/ http///aticcll//tttee//////
14

---

Published in Transactions on Machine Learning Research (12/2024)
Mario Koppen. The curse of dimensionality. 2000.
Alexander Kraskov, Harald Stögbauer, and Peter Grassberger. Estimating mutual information. Physical Review E, 6966)30 June 2004. doi: Bet.///yyyeee.e....o Publisher: American Physical Society.
URL https://link: aps.apa..... org/
doi/10. hoo//0..10.//hyseee....... 066138
Nina Kudryashova, Theoklitos Amvrosiadis, Nathalie Dupuy, Nathalie Rochefort, and Arno Onken. Paramet- ric Copula- GP model for analyzing multidimensional neuronal and behavioral relationships. PLOS Compu- tational Biology, 18(1):e100979,, Jamuary 2022. ISSN 1533- 7358. doi: 10..87//joommallllc...0000000000 URL
https://journalls plos. http://wwwppppp..tt......ao..ara
Publisher:
Public Library of Science.
Rylan Larsen, Emily Turschak, Tanya Daigle, Hongkui Zeng, Jun Zhuang, and Jack Waters. Activation of neuromodulatory axon projections in primary visual cortex during periods of locomotion and pupil dilation, December 2018. URL httpss//www. biorxiv. http://///oootet//./000/////////
Pages: 502013
Section: New Results.
Jian Ma and Zengqi Sun. Mutual Information Is Copula Entropy. Tsinghua Science Technology, 16 (1):51-54, February 2011. ISSN 1007-0214. doi: 10.1011/5100700011111000033 URL https://www.
sciencedirect. http:////.....ccccc...../c/////a
Alexander Mathis, Pranav Mamidanna, Kevin M. Cury, Taiga Abe, Venkatesh N. Murthy, Mackenzie Wey- gandt Mathis, and Matthias Bethge. DeepLabCut: markerless pose estimation of user-defined body parts with deep learning. Nature Neuroscience, 21(9):1281 1289, September 2018. ISSN 5546-1726, doi: 10. 5088//415930018-2220
URL httpps://ww, 9 Publisher: Nature Publishing Group.
https://sau.nature.. wwwp:w/.attr..oom atapa//e.ccc/rrrrrreeee/////----
Number:
Robert A. McDougal, Thomas M. Morse, Ted Carnovale, Luis Marenco, Rixin Wang, Michele Migliore, Perry L. Miller, Gordon M. Shepherd, and Michael L. Hines. Twenty years of ModelDB and beyond: building essential modeling tools for the future of neuroscience. Journal of Computational Newroscience, 42(1):1 10, February 2017. ISSN 1573-6873. doi: 10.0007/110222--11-0033-- URL https://doi.org 10. 10.1007/110227-016-00222
Lazaros Mitskopoulos, Theoklitos Amvrosiadis, and Arno Onken. Mixed vine copula flows for flexible mod- eling of neural dependencies. Fronliers in Neuroscience, 16:910122, September 2022. ISSN 1662-453X. 2022.
doi: 10.3389 fnins. 2022.91012.. URL http:://ww....tttatsss..ooorrr/p
910122/ful1
Roger B. Nelsen. An introduction to copulas. Springer scries in statistics. Springer, New York Berlin Heidelberg, 2. ed edition, 2006. ISBN 978-0-387-28659-4 978-1-4419-2109-3.
A. Onken and S. Panzeri. Mixed vine copulas as joint models of spike counts
and local field potentials. December 2016. org http:ppperr/Mxeeddeeeeeooopppppp copulas-- joint spike- -ow.oooctccooos------ooo--o--ooe counts- counts-Onken- 01150199000feeeee122222222226266
URL http::/www.ennatiicccctte www.ssataniiccchll.r.....
Janelle MP Pakan, Valerio Francioni, and Nathalie L Rochefort. Action and learning shape the activity of neuronal circuits in the visual cortex. Current Opinion in Newrobiology, 12:88-97, October 2018. ISSN 0959-4388. doi: 10. http://coob.201..0.....
011/80009448818300400
Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, Alban Desmaison, Andreas Kopf, Ed- ward Yang, Zachary DeVilo, Martin Raison, Alykhan Tejani, Sasank Chilankurthy, Benoit Steiner, Lu Fang, Junjie Bai, and Soumith Chintala. PyTorch: An Imperative Style, High-Performance Deep Learning Library. In Advances in Neural Information Processing Systems, volume 32. Cur- ran Associates, Inc., 2019. URL https://proceeddn falca2288ee//99P22ee7777777777 -taarratrtttte..ttaatttttttet.. 1tml html
15

---

http://eedenn...eerrp.....//.//e Published in Transactions on Machine Learning Research (12/2024)
Felix Pei, Joel Ye, David Zoltowski, Anqi Wu, Raced H. Chowdhury, Hansem Sohn, Joseph E. O'Doherty, Krishna V. Shenoy, Matthew T. Kaufman, Mark Churchland, Mehrdad Jazayeri, Lee E. Miller, Jonathan Pillow, II Memminng Park, Eva L. Dyer, and Chethan Pandarinath. Neural Latents Benchmark '21: Evaluating latent variable models of neural population activity, Jamuary 2022. URL http://arxiv.orrr org/ abs/2109. 04463
ar arXiv:2109.0446 04.090 (cs, q-bio].
Gregor Pirs and Erik Strumbelj. Bayesian Gaussian process factor analysis with copula for count data. Expert Systems with Applications, 197:116645, July 2022. ISSN 0957-4174. doi: 10. (0..1066).eswa.22222...... www..2022.116600.0.2 URL
https://www htttp::/www.cceeneddreeet www.sseceennddeee........oo http://wwwaaaaaaeoaooocaa//aa///
Rodrigo Quian Quiroga and Stefano Panzeri. Extracting information from neuronal populations: information
theory and decoding approaches. Nature Reviews Neuroscience, 10 185, March 2009. ISSN 1471-
0048. doi: 10.1038/ 10.1038/mm2278. nr227888. URL https: //ww: nature. com////rrrrtiiccls//mm/2778 Number: 3 Publisher:
Nature Publishing Group.
Carl Edward Rasmussen and Christopher K. I. Williams. Gaussian Processes for Machine Learning. The
MIT Press, November 2005, ISBN 078-0-262-22683- doi: 10. 10/755//mitprrss/220000.... URL https:
//direct.m mit. htt0:///.aaagaaa/////4///////eer
Sumanta Ray, Snehalika Lall, and Sanghamitra Randyopadhyay. CODC: a Copula-based model to identify differential coexpression. npj Systems Biology and Applications, 6(1):1 13, June 2020. ISSN 2056-7189, doi:
10.1038/s44 0000///ss1110-0200133--- URL
Publisher:
Nature Publishing Group.
Benjamin Rhodes, Kai Xu, and Michael U. Gutmann. Telescoping Density-Rati Estima- tion. In Advances in Neural Information Processing Systems, volume 33, pp. 4905-4916.
Curran Associates, Inc., 2020. URL Abs.raatttt..ttttt....... html
Christian P. Robert and George Casella. Monte Carlo Integration. In Christian P. Robert and George Casella (eds.), Monte Curlo Stalistical Methods, pp. 71 138. Springer, New York, NY, 1999. ISBN 978-1-4757- 3071-5. doi: 1141077998---4457--80
3. URL https:
Houman Safaai, Arno Onken, Christopher D. Harvey, and Stefano Panzeri. Information estimation using nonparametric copulas. Physical Review E, 9855)055330,, November 2018. ISSN 2470-0045, 2470-0053 doi: bio, stat].
URL http://aardt..or///a///00.....
08018 arXiv:1807.0800 [cs, math,
Eric Schulz, Maarten Speckenbrink, and Andreas Krause. A tutorial on Gaussian process regression: Mod- elling, exploring, and exploiting functions. Journal of Mathematical Psychology, 85:1 16, August 2018. com/science/
ISSN 0022-2496. doi: 10.1016/j.jmp.2011..283 atricleepp./0002249177700
URL http::://wwwws.ciinee
Gideon Schwarz. Estimating the Dimension of a Model. The Annals of Slatistics, 6(2):461-466 1978. Publisher: Institute of Mathematical
ISSN 0090-5364. URL
Statistics.
Babak Shahbaba, Bo Zhou, Shiwei Lan, Hernando Ombao, David Moorman, and Sam Bchseta. Semi- parametric Bayesian Model for Detecting Synchrony Among Multiple Neurons. Newral Computation,
(6(9))2025-205) September 2014. ISSN 0899-7667, 1530-888X. doi: 10.1162/NECO 00631. URL
https://direct..r .it. waa..eaaaancco/aaaiaa//////////0 2051//800
Claude Flwood Shannon and Warren Weaver. The mathemalical theory of communication, University of Illi-
nois Press, Urbana, 1998. ISBN 978-0-252-72546-33 URL
fy11166/8233924-bbhmm OCLC: 40716662.
Akash Srivastava, Scungwook Han, Kai Xu, Benjamin Rhodes, and Michael U. Gutmann. Estimating the Density Ratio between Distributions with High Discrepancy using Multinomial Logistic Regression, May arXiv:2305.00866 [stat].
2023. URL hTRL:tt///// hhtt:://aaxxet http://aaxi...rg///////000.0.0.0
16

---

Published in Transactions on Machine Learning Research (12/2024)
Nicholas A. Steinmetz, Cagatay Aydin, Anna Lebedeva, Michael Okun, Marius Pachitariu, Marius Bauza, Maxime Beau, Jai Bhagat, Claudia Böhm, Martijn Broux, Susu Chen, Jennifer Colonell, Richard J. Gardner, Bill Karsh, Fabian Kloosterman, Dimitar Kostadinov, Carolina Mora-Lopez, John O'Callaghan, Junchol Park, Jan Putzeys, Britton Sauerbrei, Rik J. J. van Daal, Abraham Z. Vollan, Shiwei Wang, Marleen Welkenhuysen, Zhiwen Ye, Joshua T. Dudman, Barundeb Dutta, Adam W. Hantman, Ken- neth D. Harris, Albert K. Lee, Edvard Moser, John O'Keefe, Alfonso Renart, Karel Svoboda, Michael Häusser, Sebastian Haesler, Matteo Carandini, and Timothy D. Harris. Neuropixels 2.0: A miniatur- ized high-density probe for stable, long-term brain recordings. Science 372 April 2021. doi: 10.1126 /1126//cience.abf555. URL https://www http://www.wwncceee/ee./ science. htt///////aal11.111/////e/////// abt.:/eeeeeeeeeeee..... Publisher: American Association for the Advancement of Science.
Ben Van Vliet. Abe Sklar's, March 2023. URL https:// ssrn. wtt.c//bbstraat-419955555555
Mai-Anh T. Vu, Tilay Adah, Demba Ba, György Buzsáki, David Carlson, Katherine Heller, Conor Lis- ton, Cynthia Rudin, Vikaas S. Sohal, Alik S. Widge, Helen S. Mayberg, Guillermo Sapiro, and Ka- fui Dzirasa. A Shared Vision for Machine Learning in Neuroscience. Journal of Newroscience, 38((): 1601-1607, February 2018. ISSN 0270-6474, 1529-2401. doi: 10.1523/J (5553///EUROSSL10000000
URL
https: htnpu:/// jweujnnnii....... wtt.:o//oooott////////0000 Publisher: Society for Neuroscience Section: Tech- Sights.
Sumio Watanabe. rion
PDF]
widely August 2012.
A
applicable
Bayesian information crite-
Semantic
Scholar,
URL https://www. wttp:/wwwwsaattcccooo.r..aara http:///ww.sennatcicccoooo.
A.
Wilson and H. Nickisch. Processes (KISS-GP).
Kernel Interpolation for Scalable Structured Gaus-
URL cettp::/www.semanttiie orge/ppprrr//nnnnee-ooooooorooor -et.laaal--trutturrree--eae.eeee Wilson-Nickisch 0004002acc0rrcaaaaaa0000eeeeeeee
sian
March 2015.
Byron M. Yu, John P. Cunningham, Gopal Santhanam, Stephen I. Ryu, Krishna V. Shenoy, and Mancesh
Sahani. -aussian-Process Factor Analysis for Low-Dimensional Single-TTiie Analysis of Neural Population Activity. Journal of Nearophysiology, 102(1):614 635, July 2009. ISSN 0022-3077, 1522-1598. doi: 10.
URL htap://........ll....////r//r//r
http://wwww...iooooooo..o... htapp/iollg...rrroo///........ http://....1122/......... 90941,2008
thi-Hua/ue Zhou. Why over sver-parameterizzatonn deep neural networks does not overfit? Science China Informa-
tion Sciences, 64(1):116101, Jamuary 2021. TSSN 1674-733X, 1869-1919. doi: 10.1007/s11432-D00 012114432-00--25--222 2220-288-6.
URL https://link springer www.i/nnrrr...m.m//e///oe 1011000.0000011300000 2885-6
Ariel Zylberberg. Manuel Oliva, and Mariano Sigman. Pupil Dilation: Fingerprint of Temporal Selec- tion During the " "Attentional Blink". Frontiers in Psychology, 3, August 2012. ISSN 1664-1078. doi: 10.3389 (ttps9///pggg.1122.00 2012.00316. URL
3389/fpsyg. /t1.:/pooy...211220011/1/1//
Publisher: Frontiers.
A Machine Specifications
Experiments were ran on machine utilizing 128 GB of system memory, Intel Xeon Gold 6142 processor, NVidia 2080, and NVidia 2080Ti. Plots were made a laptop with 8 GB of system memory and an Apple M2 processor. The python version utilized for all computations and plots made was 3.10.6; Elephant 1.0.0 and Gopula-GP 0.0.5 were utilized for computations, and Matplotlib 3.6.1 and Seaborn 0.10.0 were utilized for plot production.
B Information Theory Essentials
Information theory, in broad terms, is the study of formal information quantification and communication Shannon & Weaver 1998). The key measure of information in this field is entropy, defined as H(X)
17

---

Published in Transactions on Machine Learning Research (12/2024)
-Ex (log)(=)), 2 X. Essentially, entropy measures the mumber of bits of information needed to encode the total uncertainty of event X (Shannon & Weaver 1998). and as such acts as quantification of uncertainty. We also define the conditional entropy of X given the outcome of event Y as
H(XY)) Exx p(e,))), I X, V~Y, p(y)
(14)
and the matual information between X and Y as I(X: Y) = (X) - H(XYY) Mutual information thus represents the reduction in bils of uncerlainly surrounding X given the oulcome of Y. In addition, this quantity is symmetric (Shannon & Weaver 1998), and so [(X: : Y) = H(Y) H as well. In other words, mutual information describes statistical co-dependence between events X and Y. In addition, we note that mutual information is positive semi-definite (2 0), and so only ever represents a reduction in uncertainty (this is intuitive; the unconditioned outcome of an event. naturally encompasses all conditions). We extend this to the conditional mutual information of X and Y given the outcome of event Z, [(X: Z| = H( (XX,Y|Z) where (XX, I is the conditional joint entropy of X and Y given the outcome of Z, or the bits of uncertainty in the joint system of X and Y under condition z (Shannon & Weaver 1998). Finally, mutual information is casily extendable to a definition of joint mutual information via replacing singular entropy with it's joint definition (as traditional entropy itself is not the primary concern of this paper, we will skip rigorously defining these terms for brevity) (Shannon & Weaver 1998). .Jinnt mutual information and it's conditional counterpart are the primary quantities of interest this paper wish to extract the interaction information between X and Y I(X + Y) = Y-= ((XIY) -XX). representing the change in the mutual information between marginals of X when Y is learned, or the amount of uncertainty surrounding X captured in Y.
wish to estimate in
C Python Version and Hardware Utilized
Computations were made on the Edinburgh University compute server utilizing 128 GB of system memory, an Intel Xeon Gold 6142 processor, an NVidia 2080, and an NVidia 2080Ti. Plots were made on a laptop with GB of system memory and an Apple M2 processor. The python version utilized for all computations and plots made was 3.10.6; Elephant 1.0.0 and Copula-GP 0.0.5 were utilized for computations, and Matplotlib 3.6.1 and Seaborn 0.10.0 were utilized for plot production.
Table 2: Different GPLink functions used in Copula-GP for different copula. types. GPLinks are used to parameterize dependence of the particular copula in continuous variable. Courtesy of (Kudryashova et al. 2022
Copula
Domain
GPLink(f): R + dom(cy)
Independence
Gaussian
[-1,1]
Er((f/1.4)
Frank
(-o0, 00) 00, 00)
0.3 f + sign(s): (0.3 f)²
Clayton
Exp(0.3
Gumbel
(1,00
1 + Exp(0.3
D Data Retrieval and Cleaning
Visual Coding Newropizels split into recording sessions, where during each session a mouse test subject exposed to varying visual stimuli, and is publicly available via the allensdk python package (visit https://portal. brain-map. org/ and navigate to "circuits and behavior", then "Neuropixels"). Stimuli were presented in "blocks" of similar kinds of stimuli of varying duration. For this project, we examine specifically the first 100 stimuli presentations of the first "drifting gratings" block (block 2), as this way we mitigate variation in neuronal response due to differing stimulus length. In addition, within this block individual presentations are of uniform length (2 seconds of presentation, with a inter-presentation break of second), and so the block casily divisible into separate spike train trials (the GPFA implementation used is best suited for -ritialy-trial format).
18

---

Published in Transactions on Machine Learning Research (12/2024)
The session data contains the full spike recordings of all recorded units within the subject's visual cortex. As stated prior, these units possess quality metrics that correspond to how accurate and noisy unit recordings are. Of these, we utilize Signal to Noise Ratio (SNR) and Inter Spike Interval (ISI) Violation rate. SNR corresponds to the ratio of the maximum unit waveform amplitude to one standard deviation of the waveform, and acis as a metric of how noisy the unit recording is (Biau & Scornet 2016). ISI violation rate is the percentage of unit spikes that occur during what should be the corresponding neuron's 's refractory period, and serves a metric to determine whether multiple neurons and/or electronic interference are being picked up in a single unit's recordings (Czanner et al. 2015).
Session data used here was pulled from a single session (session ID 756029989). By setting lower bound for SNR and an upper bound for ISI violation rate, we are able to filter overly noisy and unreliable units' spike recordings out the data set. Utilizing a SNR lower bound of 3 and an ISI violation rate upper bound of 0.05 (5%), we found that in this session 26.0% of units meet these quality thresholds (178 out of the original total of 684). Example spike data for a single stimulus presentation of seconds is shown in figure
We also utilized pupil area recordings (cm?) to parameterize the copulas in pupil dilation. Measurements were recorded every 33ms and possess a relatively large range, with sporadic large spikes in pupil dilation (see pupil dilation curves in figure 10a). As such, we applied a rolling mean with window of 10 entries for smoothing followed by a robust normalization procedure Roobst(X):
X- (x))
Robust Qx(X) - Qx
(15)
where Qx)), Qx( are the 1st and 3rd quartile values of X. Doing so, we reduce the impact of outliers in the data. As the input for parameterizing values for fitting C-vine via Copula-GP must be the interval [0, 1}, we map onto via an additional min-max normalization
min(X)
MinMax(XX= = max (X) min(X)
(16)
We can observe the difference between the raw and processed data in figure 10 See figure 10c for the distribution this preprocessing regime created, which appears roughly normal with few outliers.
The technical white-paper for the dataset is also available at:
divio ddvoommedia..nt///o1o__.011/
75/ 71/00//70000000----030--0--66220 liaiiiieeesl__ssslleeeeeaaeee_le visual_coding. - __white white_ppperr v10. pdf
E
Copula-GP Implementation
The vine construction implemented in Copula-GP used copula building blocks from five distinct families: independence, Gumbel, Gaussian, Frank, and Clayton copulas, with the independence copula preferred for independent variables (see figure Gb. The factorization of the vine is that of the C-vine described by equation & with accommodations made for EP-parameterization. To fully capture the tail dependencies and negative correlations in relationships between marginal distributions, mised copulas were utilized. In the core paper, these are defined as
K Cmiizdd X(V) Σ oy (XC,(1
(17)
where K is the number of elements, dj is the concentration of the jth copula (cj), and 02 is the jth copula's 's parameter GPLink. The GPLink for each copula is determined by it's copula variant, as shown in figure 2 with 0 being defined by GPLink(J), where f is sampled from ~ N(u, (x(X,)) (the choice of GPLink depends on the kind of copula; see 2.. GP is also utilized to parameterize the concentrations Φ3, which are defined as
و= (1 (1-11-1 (11-)) (g) II bm Im = & (Im + -1
MM- - m -----T--- -)
IM = 0,
(18)
M m
19

---

Published in Transactions on Machine Learning Research (12/2024)
where @ is the CDF of standard normal distribution and fm is sampled from fm Niim: K (Y) This gives us 2M -1 sets of hyper parameters to estimate, [X] M kernel hyperparameters for each GPLink 0 and (A}M-1 kernel hyperparameters for each concentration function 6, estimated via the methods described in section
E.1 Copula-GP Model Selection and Parameter Estimation
As stated in the body of this paper, the parameters for the distributions used in GP must be estimated. Tn the Copula-GP framework, this is accomplished through the use of stochastic variational inference (SVI), with SVI being scaled to high dimensions by mcans of Kernel Interpolation for Scalable Structured Gaussian Process (KISS-GP) (Wilson & Nickisch, 2015). These methods were specifically used for efficient implemen-
tation in aforementioned python libraries Pyforch and GPyTorcch (Paszke et al. 2019 Gardner et al. 2018 For model hyper-parametee parameter selection, Watanabe- Akaike information criterion (WAIC) was used, a metric which aims to maximise the Akaike information criterion (AIC) by means of a Bayesian-approach (we more rigorously define ATC in section 7) Given f1, f2 and the true distribution 9, WAIC examines the difference
in log- likelihhoodd i.e. [g: fi) [g: f2) = (Eolo (i()) log f2(X)), and selects the model with the greater mutual information (Watanabe 2012). The space of models is too large to find the optimal model when considering the number of combinations of copulas shown in 6 and as such the core paper implements a greody algorithm of minimising WAIC which can be used with all copula types and a heuristic algorithm specifically tuned for certain combinations of copula types. For this project, we utilize the heuristic approach of copula selection.
F GPFA Specifics
F.1 Dimensionality Selection
Before we can utilize GPFA, we first must isolate which number of dimensions n to reduce down to. For the in vivo dataset we chose n via investigating the log likelihood of GPFA fits extracted from a 3-fold cross validation, from target dimensionality n to 50. We then plotted the log likelihoods and saw both where the clbow, or the point of maximum concave slope curvature representing a "good enough" dimension to reduce down to (Antunes al., 2018) of the log likelihood curve was. While there are computational methods for elbow selection (Antunes et al. 2018 to be selected ourselves visually.
investigated few enough points to allow for the elbow
F.2 Additional Post-GPFA Interim Processing Required
Copula-GP's C-vine framework fit on single-trial continuous, however Elephant's GPFA implementation produces trajectory data that is split into trials. A solution this would be to concatenate trials trajectory wise, however per-trial drift in trajectory means result in weak Copula-GP fit performance if these lead to large jump discontinuities and thus loss of smoothness in the data; the kernel for the GP-link functions used in parameterization will require more restrictions as it encodes the smoothness of the data (among other things) (Schulz et al. 2018 Yu et al.) 2009). For the in vivo data, we found in figure 7 that drift occurs in the second inter-stimulus break in the average trial. As such, we crop this period out of cach trial (50 points), and concatenate trials logether trajectory-wise. While this by far not perfect solution, it allows the data to remain roughly smooth at the cost of continuity time and residual (small) jump discontinuities, as well isolating the data only when stimulus presentations are occurring. See figure for example of tial-to-trial discontinuities created by this interim step.
g
F.3 Python Implementation of GPFA Used
The implementation of GPFA utilized is sourced from the Elephant (Electrophysiology Analysis Toolkit) python library, with the package's 1.0.0 release (used for the contents of this paper) being published November 10, 2023 (Denker & Kern 2023). The Elephant package was specifically designed for use on neuronal data, motivated by push to release standardized python package for in computation neuroscience. The
20

---

Published in Transactions on Machine Learning Research (12/2024)
GPFA module has specifically use in recent papers Pei et al. 2022 Bagi et al. 2022 and since release has become quite popular. The GPFA module receives some dimension to reduce down to and time-bucket size m (we utilize 10ms and 20ms buckets) to instantiate a GPFA python object. This object can then be fit on number of spike train recordings of uniform temporal length given the start and end time of cach recording, summing the number of spikes for each time bucket and utilizing the bucketed spike counts as the observation matrix equation equation 13 As each projection corresponds to the expectation of the latent trajectories [(O (), these parameters estimated in the Elephant implementation via expectation maximization Yu et al. 2009).
Algorithm 1 Implemented bivariate mixture copula. bagging algorithm. Here, mixture copulas consist of mixed copula. variants variants, their corresponding dependency parameter values θ (for each of the m inputs), and their corresponding mixing parameters már. We assume weights have been defined over all m input points.
1: Input
mixture: A list of mixture copulas aggregate. weights: Corresponding mixture weights.
4: Output
5: míature Variant List, Mix Parameters, θ): A new mixture copula of variants Variant List with correspouding mixture parameters Mix Parameters and dependency parameters weights)
6: procedure BAGCOPULAS
7:
N
b. Unique variant counter.
8: indexes +- Emply dictionary
b Indexing dictionary.
9: Variant Counts + Empty dictionary
> Counting dictionary.
10: Total Variant Weight + dict()
b. Weighting dictionary.
11:
for each i-th mixture copula in mixlures do for each n-th variant in the mixture copula do if the variant has not yet been encountered then Let the variant act as a key for N from indezes Let the N-th index of Vorianl Counls get Let the N-th index of Total Variant Weight get Increment N by 1
12:
13:
14:
15:
16:
17:
18:
idr +- corresponding entry for the variant in inderes Increment the idx-th entry of Variant Counts by
19:
20:
Let (i, n) act as a key for idæ from indexes
21:
Add the the i-th entry of weights to the idx-th entry of Total Variant Weight
22:
Variant List +- an emply list Mir Parameters + ((Nx))
23:
> Mixture parameters.
24:
θ - 0
Dependency parameters.
25:
for each i-th mixture copula in miatures do for each n-th variant in the mixture copula do
26:
27:
idz - corresponding entry for (i,n) n) in indezes Set the idæ-th entry of Variant List the variant weight t- the i-th entry of weights
28:
29:
30:
+ the n-th mixture parameter of the mixture Add x weight to the idr-th row of Mir Parameters
31:
32:
+ the n-th dependency parameter the mixture Total Weight <- the idr-th entry of Total Variant Weights Add weight/ Total Weight to the idx-th row +
33:
34:
35:
return A new mixture copula mixture Variant List, Mir Parameters, +)
21

---

Published in Transactions on Machine Learning Research (12/2024)
G
Validation Test Reproducability
All validation tests were done using variations of arguments for the validate validate_bagging py python file.
Test 1 For test 1, we utilized:
Random Seed: 859448723
Num. Estimators: 4
Maximum Copula Elements: 5
Dimensions:
Shuffling of Data: Yes
Input Type: Random
Test 2 For test 2, we utilized:
Random Scod: 859443
Num. Estimators: 4
Maximum Copula Flements: 3
Dimensions:
Shuffling of Data: Yes
Input Type: Random
Test 3 For test 3, we utilized:
Random Seed: 859443
Num. Estimators: 4
Maximum Copula Elements: 3
Dimensions:
Shuffling of Data: No
Input. Type: Linear
22

---

Published in Transactions on Machine Learning Research (12/2024)
H Additional Figures
Gaussian E [--1111
Frank 0 E
Clayton è [0,
Gumbel e [1. m))
Clayton 90" Gumbel
1.0
1.0
1.0
1.0
1.0
0.8
0.8
0.8
0.8
0.8
0.6
0.6
0.6
0.6
0.2
0.2
0.2
0.2
0.2
6
0.8
0
10.0
0
5.0
5.0
0
0.0
0.0
0.0
0.0
0.0
1.0
0.5
1.0
0.5
1.0
0.5
y1
yi
Mixture of copulas copula
Rotation-equivariant copulas
Copulas that can be rotated by 90
Figure 6: Various kinds of copulas utilized in Copula-GP. Note the difference in tail distribution represen- fation, as well as how combining copula variants make a new mixed copula. Original source (Kudryashova
al. 2022 with permission.
Trajectories (Min-Max Scaled)
Average Trajectories (Min-Max Scaled) per Stimulus Trial
over time for 2 continuous trials
0.6
0.8
0.5
0.6
0.4
0.4
0.3
0.2
0.2
Start of stimulus presentation
0.1
End of stimulus presentation
0.0
0.0
0.5
1.0
1's
2.0
2.5
3.0
Time (s) since Stimulus Onset
Time (s) Stimulus Onset
(a)
(b)
Figure Trajectories extracted via GPFA (dimension n 13) representing dynamics driving recorded neuronal activity. Figure 7a are two trials' trajectory data as functions time that have been concatenated, giving the appearance continuity. Figure 7b is the mean trial. Note the clear vertical drift in mean and single trial trajectories post-stimulus presentation stop. Data is scaled to the range [0.01,0.99] to match G input range of (0,1). Lines transluscent for visibility purposes.
23

---

Published in Transactions on Machine Learning Research (12/2024)
Time Discontinuity in Trajectories Post Cropping
Break point
0.50
0.45
0.40
0.30
0.25
0.20
60
80
100
120
140
Figure 8: Discontinuities present in trajectory data post-cropping and -concatenation (3 of 13 trajectories shown). Such discontinuities can negatively effect GP performance when estimating GP-link functions during the Copula-GP fit process.
Log Lik. for Num. Dimensions n
Highest Log Lik.
Elbow
Log-Likelihood
-5000
-10000
-15000
10
13
37
40
43
46
49
Dimensionality
Figure 9: Log-likelihood curve for GPFA fit on trials as function of dimension n. Elbow and max log-
likelihood found marked. Note that only dimension up to and including n 50 were tested, and entirely possible likelhold increase further in dimensions.
24

---

Published in Transactions on Machine Learning Research (12/2024)
Pupil Area over All Trials (600s total) Raw
Pupil Area over All Trials (600s total) Normalized
1.00
Pupil Area (cm²)
Pupil Area (Normalized)
10-2
0.75
0.50
0.25
10-3
0.00
200
400
600
200
400
600
Seconds
Seconds
(a)
(b)
Distribution of Pupil Values over Full Trials
Distribution of Pupil Values over 15 Continuous Trials
Mean: 0.5478
101
Q1: 0.5340
101
Q3:
Count
Count
10
10
Mean: 0,5480
Q1: 0.5353
Q3: 0.5636
0
0.0
0.2
0.4
0.6
0.8
1.0
0.50
0.55
0.60
0.65
Pupil Dilation (Normalized)
Pupil Dilation (Normalized)
(d)
Figure 10: Figures 10a and 10b represent pupil dilation function of time over all trials (duration 600s). raw 110a pupil area in em2 scaled logistically) and processed 10b pupil area smoothed and normalized). Processing consisted of rolling mean (window of 10 entries), robust scaling, and a min-max scaling, Note the removal of strong outliers through processing, as well as the large range necessitating a logistic scale in values present in the raw recordings. Figure 10c represents the distribution over all pupil dilation data utilized (100000 continuous samples / 100 trials). Figure 10d represents the distribution over 1500 continuous samples (15 trials). Note the use of logistic scale in figures 10al
10c
and 10d
25

---

Published in Transactions on Machine Learning Research (12/2024)
spike raster for stimulus presentation 3997
unit
12
2183.25 2183.50 2183.75 2184.00 2184.25 2184.50 2184.75 2185.00 2185.25
time (s)
Figure 11: Spike raster of a single drifting gratings stimulus presentation (duration of 2s), from the session utilized model validation (session ID 756029989). Only those spikes with SNR lower bound and ISI violation rate upper bound of 0.05 are shown. Note correlation in spike events present the raster plot, well as heightened neuronal spike-rate at the beginning of stimulus presentation.
26


Reviews:
[{"review_id": "4LFMfi5rxL", "reviewer": "Reviewer_172m", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Strengths:\n* The paper is well written.\n* The paper includes a simple introduction to a Coupla-GP and clear explanation of bagging.\n* The paper is a comprehensive application of a bagged Copula-GP (and baselines) for pupil dilation dependence on visual cortex trajectories. I especially appreciate the naive average.\n\nWeaknesses:\n* \"GP is in essence a regression attempting to estimate a relationship between variables by means of Bayesian inference.\" This statement is significantly over simplified. A GP can be used for regression but is not a regression. I would also suggest that you cite the GP textbook here (Rasmussen & Williams, 2006).\n* There is little discussion of why there is a need for an estimate approach to Copula-GP. \"Computing this directly via integration over the estimated copula density is computationally expensive due to nested integrals.\" A runtimes comparison for both approaches should be included (at least for the toy dataset) with a more detail discussion of infeasibility for the *Neuropixels* dataset.\n* I would like to see standard deviations for results in Table 1. It is unclear how important a difference in the thousandth place is without them.\n* Input and label details for the *Neuropixels* dataset could be clearer.\n\nCarl Edward Rasmussen and Christopher KI Williams. *Gaussian Processes for Machine Learning*. The MIT Press, 2006.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "Yes", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "This paper applies existing techniques that shed light on the strengths and weaknesses of Copula Gaussian process (GP) based methods for pupil dilation dependence on visual cortex trajectories. They include results on toy dataset that show dynamically aggregated copulas via Bayesian Information Criterion (BIC) maintain competitive performance while reducing computational costs.", "claims_and_evidence_yn": "Yes", "claims_explanation": "", "audience_yn": "Yes", "audience_explanation": "", "requested_changes": "* Paper readability would benefit from using `\\citet` and `\\citep` latex commands.\n* In Hu & O’Hagan (2021), $N$ is used to represent the $|X|$ whereas in the current work $N$ is the number of Monte Carlo (MC) samples. Consider using $S$ to represent the number of MC samples.\n* Axes for Figures 1, 2, and 3 would be helpful.\n* I would like to see the truth or mean plotted for Figures 4 and 5.", "broader_impact_concerns": "N/A", "additional_comments": "", "extra_scores": {}, "extra_text": {}}, {"review_id": "fAuEDR4eVg", "reviewer": "Reviewer_CXvm", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "## Strengths\n- Motivation and Background: the paper discusses the challenges in estimating the statistical dependence in neuronal data, where it is typical to deal with high-dimensional observations. \n-  The empirical validations are sound and help demonstrate the effectiveness of the proposed method.\n\n## Weaknesses\n- The paper needs a major revision on the writing: technical background sections (Section 4-5-6) need to have more rigorous formalism -- for the moment they are too verbose; citations get duplicated in different places, and almost all of the citations are done in the wrong way (need to use \\citep); the authors should break large paragraphs into smaller ones. Overall, in my opinion, the writing and formalism, in the current version, do not match the usual standard of a Machine Learning venue.\n-  The novelty is quite limited. I believe all of the core methodologies in the papers are taken from previous works. For example, the proof of mutual information equals copula entropy is the same as in the original paper that introduced them (Ma and Sun 2011). The authors, however, did indeed state that their submission is a \"validation of combining Copula-GP with GPFA\", adding the bagging as the ensemble.\n- I have a concern about the high dimensionality of the realistic dataset benchmark in Sections 9 and 11: in the main text, the original dimensions of the problem were not mentioned.\n- The runtime of the proposed methods is not mentioned in the main text. The method might be effective, but if the runtime is too prohibitive, its practicality will remain in question. \n- Nitpick: algorithm 1 is written in a pseudo-python style. It should be more generally written for readers who use other programming languages in data science/Machine Learning to understand (such as R or Julia).", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "Yes", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "This paper introduced a framework of ensemble multiple copula estimators, with the goal of estimating mutual information/information interaction between groups of brain's neurons conditional to behavioral variable. In particular, the authors proposed to combine bagging methods with dimension reduction using Gaussian Process Factor Analysis (GPFA) and Copula-GP, and show with numerical validations on simulated and realistic (neuronal data of visual cortex on pupil dilation) datasets that this combination of frameworks help to improve the performance of the mutual information estimation task.", "claims_and_evidence_yn": "Yes", "claims_explanation": "", "audience_yn": "No", "audience_explanation": "", "requested_changes": "-  I suggest the authors to move many parts of the formalism in the Appendix to background sections 4-5-6.\n- Math typo in the first paragraph of Section 6: mutual information $I(X, Y) = H_C(X) - H_C(X \\mid Y)$, not the other way around. \n- I am aware of two of the works that can be used to estimate Mutual Information that can scale to high-dimensional data: TRE (Rhodes et al. 2020) and DRE (Choi et al. 2022). Could the authors perform the mutual information estimation task with high-dimensional Gaussians in Section 6.2 of Choi et al. (2022) with their framework included?\n- Could the authors provide the runtime of the methods, separated by each step (GPFA, fitting copula, ensembling)?\n\nRhodes, Benjamin, Kai Xu, and Michael U. Gutmann. \"Telescoping density-ratio estimation.\" Advances in neural information processing systems 33 (2020): 4905-4916.\n\nChoi, K., Meng, C., Song, Y., & Ermon, S. (2022, May). Density ratio estimation via infinitesimal classification. In International Conference on Artificial Intelligence and Statistics (pp. 2552-2573). PMLR.", "broader_impact_concerns": "N/a", "additional_comments": "", "extra_scores": {}, "extra_text": {}}, {"review_id": "OgKQEVT6jm", "reviewer": "Reviewer_paFS", "summary": "", "questions": "", "limitations": "", "rating": null, "confidence": null, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "**Strenghts**\n- The paper presents an extensive experimental section in which the influence of the proposed additions to Copula-GP are evaluated in contrast with the lack of them.\n- The paper provides sufficient background information on what a Copula is, what Copula-GP is to make the paper easier to follow to readers that are not familiar with the topic.\n\n**Weaknesses**\n- The presented work looks more like a set of proposed additions over Copula-GP to improve its performance rather than a simple clear idea. \n- The dimensionality of the data is not very high. \n- The paper's structure has too many sections, some of them being very short. The paper would benefit from a re-structuring, by adding more subsection and reducing the number of sections.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "Yes", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "paper_topic_and_main_contributions": "", "excitement": null, "reproducibility": null, "ethical_concerns": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "summary_of_contributions": "The presented work introduces two additions on top of Copula-GP to improve the estimation of the dependency between a set of variables and tested on neuronal data.\n\nOne of the proposed additions is combining GPFA with Copula-GP. GPFA is a dimensionality reduction strategy that helps the Copula-GP to perform in a lower dimensional space.\n\nThe other proposed addition is using bagging in combination with Copula-GP. Bagging is a bootstraping strategy in which several 'weak' models are trained and the prediction of all models is combined to provide a stronger prediction. The authors propose learning multiple copulas and provide the estimated copula as a the mean of the mixture.", "claims_and_evidence_yn": "Yes", "claims_explanation": "", "audience_yn": "Yes", "audience_explanation": "", "requested_changes": "- Improve Structure of the paper: Adapt the paper to have less larger sections, rather than a lot of short sections.", "broader_impact_concerns": "None that I am aware.", "additional_comments": "", "extra_scores": {}, "extra_text": {}}]

Author rebuttal:

Meta-review: The paper is found to be sound, with a relatively narrow scope.

* Please share openly the code.
* it's copula variante -> its copula variant
* "As such, in addition to a naive average, we may choose to weigh models based off various Bayesian criterion." -> "As such, in addition to a naive average, we may choose to weigh models based on various Bayesian criterion."
* "it’s parameterization" -> "its parameterization"
* "catch dependencies in distributions that the baseline perceives as independence." -> "catch dependencies in distributions that the baseline perceives as independent."
Decision: Accept (poster)