Bayesian inference for low rank spatiotemporal neural receptive fields
The receptive field (RF) of a sensory neuron describes how the neuron integrates sensory stimuli over time and space. In typical experiments with naturalistic or flickering spatiotemporal stimuli, RFs are very high-dimensional, due to the large number of coefficients needed to specify an integration profile across time and space. Estimating these coefficients from small amounts of data poses a variety of challenging statistical and computational problems. Here we address these challenges by developing Bayesian reduced rank regression methods for RF estimation. This corresponds to modeling the RF as a sum of space-time separable (i.e., rank-1) filters. This approach substantially reduces the number of parameters needed to specify the RF, from 1K-10K down to mere 100s in the examples we consider, and confers substantial benefits in statistical power and computational efficiency. We introduce a novel prior over low-rank RFs using the restriction of a matrix normal prior to the manifold of low-rank matrices, and use "localized" row and column covariances to obtain sparse, smooth, localized estimates of the spatial and temporal RF components. We develop two methods for inference in the resulting hierarchical model: (1) a fully Bayesian method using blocked-Gibbs sampling; and (2) a fast, approximate method that employs alternating ascent of conditional marginal likelihoods. We develop these methods for Gaussian and Poisson noise models, and show that low-rank estimates substantially outperform full rank estimates using neural data from retina and V1.
In this paper, the authors develop a sophisticated set 
      of Bayesian techniques to infer low-rank neural receptive fields (RFs). 
      This is a problem of interest to many modellers and experimenters within 
      sensory neuroscience, as identifying the transformations from stimulus 
      spaces to output spike trains is both extremely useful for experimental 
      protocols, and very difficult or data-intensive due to the very high 
      dimensionality of the stimulus space. One option has been to exploit 
      potential symmetries in stimulus space, and build low-rank approximations 
      to RFs based on sums of space-time (e.g. Pillow et al, 2008) or 
      frequency-time (e.g. Ahrens et al 2008) separable kernels. While these 
      methods are useful, they have, to date, received less attention from 
      model-development front. This has limited their widespread use. 
      Here, the authors extend a set of Bayesian inference techniques to 
      low-rank RF estimation. This consists of adding hierarchical priors to the 
      components of the reduced-rank RFs (using the ALD framework of Park & 
      Pillow, 2012), providing full-Bayesian and fast approximate inference 
      algorithms for the hierarchical model under a Gaussian likelihood, and 
      fast approximate inference method for the popular linear-nonlinear-Poisson 
      (LNP) likelihood using a Laplace approximation. The authors then 
      demonstrate the utility (and speed) of these methods by fitting low-rank 
      RFs to simulated and real neural data from retina and V1. Overall, this 
      work is successful and pleasing, generally clear (constrained by the 
      required brevity of the format), original, and likely to yield promising 
      future results. My only constructive comments on the paper would 
      be to clear up a few confusions and add some minor detail: - Lines 
      133-141. My understanding of this is that there are separate 
      hyperparameters {rho, ..., Phi_f} for C_t and C_x, but this is not clear 
      within the text. The current version is written to describe ALD for 
      full-rank RFs, and should be amended for the low-rank cases. In the 
      sentence "In the ALD prior..." (line 135), the localisation of C_x should 
      be in space and spatial frequency, while that for C_t should be in time 
      and temporal frequency. Likewise, the description of mu_s and mu_f as 
      being centred in space-time and frequency should change (and line 140 
      too). Also, mu_s and mu_f as length D vectors, and Phi_s and Phi_f as DxD 
      matrices, should be D_x / D_t and D_x by D_x / D_t by D_t respectively, 
      with footnote 1 changed to reflect this. - Figs 1-2: space/time 
      labelling of axes, as in Fig. 3, would be helpful. - The results 
      in Figs 3-4 assume a Gaussian likelihood, yes? I am guessing this due to 
      the use of the Gibbs sampler. This should be stated. - Are the 
      gains in model performance shown in Figs 3-4 the result of using the 
      low-rank approximation, or the result of using the ALD prior? It would be 
      helpful to quantify how much is being gained or lost by making the 
      low-rank assumption (especially since a major contribution of this paper 
      is to extend ALD to the low-rank case). To this end, a comparison against 
      a full-rank ALD estimate would be of value: just demonstrating an 
      improvement over STA is a straw man. From Fig 5B, my guess would be that 
      low-rank ALD would perform the same as full-rank ALD for the Gaussian 
      likelihood used in Figs 3-4. Notwithstanding this, the order-of-magnitude 
      improvement in compute time for low-rank over full-rank remains useful. 
      - As always, the experimenters who would find this work most 
      useful would also likely find it near impossible to digest. The net 
      utility of this work depends on code being made 
      available.