Dynamic functional connectivity among neuronal population during modulation of extra-classical receptive field in primary visual cortex

Brain Research Bulletin 117 (2015) 45–53

Contents lists available at ScienceDirect

Brain Research Bulletin journal homepage: www.elsevier.com/locate/brainresbull

Research report

Dynamic functional connectivity among neuronal population during modulation of extra-classical receptive field in primary visual cortex Xiaoke Niu a , Li Shi a,b,∗ , Hong Wan a,∗∗ , Zhizhong Wang a , Zhigang Shang a , Zhihui Li a a b

School of Electrical Engineering, Zhengzhou University, Zhengzhou 450001, China Department of Automation, Tsinghua University, Beijing 100000, China

a r t i c l e

i n f o

Article history: Received 1 December 2014 Received in revised form 3 July 2015 Accepted 8 July 2015 Available online 17 July 2015 Keywords: Extra-classical receptive field Surround modulation Functional connectivity Primary visual cortex Granger causality

a b s t r a c t The neuronal activity evoked by stimuli confined in a receptive field can be modulated by surround stimuli of the extra-classical receptive field (eCRF). The surrounding modulation, hypothesized to be the basis of visual feature integration and figure-ground segregation, has drawn much attention in the field of neuroscience and engineering. However, most studies focused on surround modulation of individual neuronal response. In this study, we analyzed surround modulation of the population response recorded from rat primary visual cortex, and further investigated dynamic functional connectivity modulated by the surrounding stimuli. The functional connectivity was estimated using Granger causality (GC) and then determined by thresholding the p-matrix with different significance ˛ values. Four scalar indexes were calculated to describe the functional connectivity of neuronal population: averaged connection strength (mGC), connection density (D), clustering coefficient (C) and path length (L). The statistical results from 5 rats showed that these network characteristics were dynamically changed during modulation of surrounding stimuli, which suggested that the neuronal population may connect in a dynamic way during modulation of eCRF. We further guessed that the neurons may happened to be organized in a more efficient way underlying surrounding modulation conditions, which helps to process larger images efficiently with the same number of neurons. This study provided new insights for a better understanding of the underlying neural mechanisms responsible for surround modulation. © 2015 Elsevier Inc. All rights reserved.

1. Introduction The stimuli surrounding the classical receptive field (CRF), cannot elicit spikes when stimulating alone, but can modulate the visual response to CRF stimuli (Barlow, 1953; Hubel and Wiesel, 1965; Cavanaugh et al., 2002a). This phenomenon refers to the so called “surround modulation,” which is hypothesized to be the basis of visual feature integration and figure-ground segregation (Albright and Stoner, 2002; Seriès et al., 2003). Furthermore, understanding the underlying neural mechanisms is important for developing a theoretical model of early signal integration and neural encoding of visual features in the primary visual cortex (V1).

Abbreviations: RF, receptive field; CRF, classical receptive field; eCRF, extraclassical receptive field; GC, Granger causality; V1, primary visual cortex; LFP, local field potential. ∗ Corresponding author at: School of Electrical Engineering, Zhengzhou University, Zhengzhou 450001, China. Fax: +86 0371 67781407. ∗∗ Corresponding author. Fax: +86 0371 67781407. E-mail addresses: [email protected] (L. Shi), [email protected] (H. Wan). http://dx.doi.org/10.1016/j.brainresbull.2015.07.003 0361-9230/© 2015 Elsevier Inc. All rights reserved.

Basic studies on surround modulation have shown that the surround stimuli modulated both the spiking activity of individual neurons and slow components (<250 Hz) of the extracellular potentials (local field potential, LFP) (Zhang and Li, 2013), the latter of which generally reflected synaptic inputs. Further studies have found that the surround modulation, mostly suppressive, was selective to visual features, such as orientation, spatial frequency, and stimulus size (Cavanaugh et al., 2002b; Webb et al., 2005), and that surround suppression for increasing stimulus size (“size tuning”) had preferred receptive field sizes that depended on contrast (Henry et al., 2013). It was also reported that the surround suppression would sharpen the orientation tuning (Okamoto et al., 2009) or enhance the orientation selectivity (Chen et al., 2005; Xing et al., 2005) of a single neuron in V1. Thus far, the majority of studies on surrounding modulation have concentrated on exploring the modulation of individual responses, such as the spiking activity at the single neuron level (Osaki et al., 2011; Henry et al., 2013) or LFP at individual recording locations (Gieselmann and Thiele, 2008; Zhang and Li, 2013), under the surrounding stimuli. There were also many other attempts to explain the underlying neural mechanisms responsible for surround modulation and those

46

X. Niu et al. / Brain Research Bulletin 117 (2015) 45–53

Fig. 1. Examples showing implantation of a 4 × 4 microwire electrode array in rat V1, and a histological slice used to verify electrode tracks in V1.

studies mostly hypothesized that the modulation of extra-classical receptive fields (eCRF) in V1 was a product of long-range horizontal connections (<5 mm) (DeAngelis et al., 1994; Dragoi and Sur, 2000; Stettler et al., 2002; Reynaud et al., 2012), and/or feedback from extrastriate areas (Angelucci et al., 2002; Sceniak et al., 2002; Bair et al., 2003; Wang et al., 2010; Nassi et al., 2013). There were also differing reports that response modulations by static texture surround in V1 were irrelevant to the feedback connections from V2 (Hupé et al., 2001) and might have resulted from local, shortrange (<0.5 mm) lateral connections within V1 (Wielaard and Sajda, 2005). The mechanisms remain unknown. Nevertheless, it could be speculated that the surround suppression phenomenon involved interactive activities of the neuronal population, and was associated with the interaction activities between individual neurons or between neuronal groups, locally and/or globally. In our study, we focused on investigating the surround modulation of functional connectivity among the neuronal population from V1, recorded with microelectrode arrays, without regard to the spatial scale. Considering the horizontal connections were frequently observed in layer 2/3 (Burkhalter, 1989; Stettler et al., 2002; Medini, 2011), and were much less affected by long-range axonal targets from higher areas (Brown and Hestrin, 2009), we mainly concentrated on analyzing the functional connectivity in layer 2/3 of V1. With the rapid development of computational techniques for estimating the functional connectivity between or among neurons, powerful approaches have emerged (Yook et al., 2013). The commonly used statistical/computational methods were Granger causality (GC) (Achard et al., 2006; Guo et al., 2008; Ge et al., 2009; Ouyang et al., 2014), dynamic Bayesian inference (Lee et al., 2006; Neapolitan, 2009; Ide et al., 2014), the maximum likelihood model (Okatan et al., 2005), and the Ising model (Yu et al., 2008), of which GC analysis has become an increasingly popular method for identifying the causal connectivity in neural time series data (Barnett and Seth, 2011). Recently, it has also been used to reconstruct the underlying anatomical connectivity in conductance-based integrate-and-fire neuronal networks (Zhou et al., 2014). In this study, we analyzed the modulation of averaged population responses, recorded from layer 2/3 of rat V1, against the size of the visual stimuli. The CRF and eCRF for the population were further approximated using the same method applied on individual neurons. The functional connectivity among the neuronal population was computed with the GC method (Barnett and Seth, 2014), under CRF alone and CRF combined with eCRF stimuli conditions. Four characteristics were measured to describe the dynamics of functional connectivity among neuronal population: connection strength (mGC), connection density (D), clustering coefficient (C) and path length (L). The statistical results proved that the surround stimuli would dramatically changed the functional connectivity among the neuronal population from V1.

2. Material and methods 2.1. Animal preparations Our data were recorded from five Long Evans (LE) rats (2–4 months of age, 200–300 g) supported by the Animal Center of Zhengzhou University. All experimental procedures were in accordance with the guidelines of the National Institutes of Health and were approved by the Animal Care and Use Committee of Zhengzhou University. Before surgery, each rat was anesthetized with urethane (1.3–1.5 g/kg body weight) and restrained in a stereotaxic apparatus (David Kopf Instruments, Tuhunga, CA, USA). The rectal temperature was monitored and body temperature was maintained at 37.5 ◦ C with a heating pad. The right eye was fixed with a metal ring to prevent eye movement, the pupils were dilated with topical application of 1% atropine sulfate, and the nictitating membranes were retracted with 5% phenylephrine. A craniotomy (diameter >2 mm) was performed over the left V1 with the core at 3 mm lateral and 7.5 mm posterior to the bregma (Fig. 1a). The microelectrode array (MEA) was implanted into the area (marked with a black square in Fig. 1a) after being dipped in 3% DiI solution (1,1 -dioctadecyl- 3,3,3 ,3 tetramethylin-docarbocyanine perchlorate; invitrogen, Carlsbad, CA, USA) prepared with dimethylformamide (DMF) to label the track of the electrode insertions (Zhang and Zhang, 2010). The MEA consisted of 16 polyimide-insulated platinum/iridium microwires (Clunbury Scientific, Bloomfield Hills, MI, USA) that were arranged in four rows with four wires in each row (electrode diameter = 50 ␮m; electrode spacing = 350 ␮m; row spacing = 350 ␮m; impedance = 20–50 k). Then the array was lowered approximately 200–600 ␮m below the V1 surface (Fig. 1b), corresponding to layer 2/3 of V1, where horizontal connections could be frequently observed. Finally, a silver wire from the array was connected for grounding to bone screws inserted into the skull surrounding the surgical openings. The recording sites were finally examined with histological slices, under a fluorescence microscope (Fig. 1c, see Section. 2.7 for details). 2.2. Visual stimuli Visual stimuli were generated from a PC running a Matlab toolbox (Psychtoolbox; MathWorks, Natick, MA, USA), and displayed on an LCD monitor (250 × 400 mm, refresh rate, 60 Hz) positioned 15 cm in front of the rat’s right eye. The center of the receptive fields for individual neurons were determined with sparse noise stimuli, with a single bright square (6.6 × 6.6◦ ) flashing on a black background at each of 12 × 12 positions in a pseudo random sequence (20 flashes/position). The spatiotemporal structure of RF was estimated with reverse correlation (Jones and Palmer, 1987). The spatial profile at each temporal delay was approximated with

X. Niu et al. / Brain Research Bulletin 117 (2015) 45–53

47

Fig. 2. The spatial distribution of classical receptive fields of each neuron group recorded from each of five rats. The red dot in each figure indicates the center of the mosaic area for the neuron group. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

a two-dimensional (2-D) Gabor function (Sun and Dan, 2009; Liu and Yao, 2014):





−((x−x0 )cos+(y−y0 )sin)2 /x2 −((y−y0 )cos−(x−x0 )sin)2 /y2

Rs x, y,  = We

(1)

where W is amplitude from the RF response, x0 and y0 refer to the locations of the RF centers,  x and  y determine the width and length of the RF, and  is the orientation. Surround modulations upon population responses were tested under different sizes of drifting sinusoidal gratings (time frequency of 2 Hz, spatial frequency 0.1 cpd, and contrast of 100%) of 12 different orientations (0–330◦ , 30◦ in steps) in circular and annular apertures on a gray background with the same mean luminance as the gratings. The center of the stimuli was determined by the mosaic area of the RFs of individual neurons (Fig. 2). Each orientation was randomly presented as the size increased. Each stimulus was repeated for 10–15 times, 2 s in duration, with a 1 s blank window in between, to avoid adaptation. Functional connectivity among networks was estimated under the same sets of drifting sinusoidal gratings, but of two specified sizes, which were determined by the size tuning curve of the averaged population response (see Section 3.1 for details). 2.3. Electrophysiological recordings The signals were collected with a Cerebus system (Blackrock® Microsystems, Salt Lake City, UT, USA) and amplified 4000×. The spikes were extracted by band-pass filtering (Second-order Butterworth) raw signals between 250 Hz and 5 kHz with a sampling rate of 30 kHz. The LFP was extracted by band-pass filtering (Four-order Butterworth) raw data between 1 Hz and 250 Hz, and sampled at 2 kHz.

chronux,org/, in this procedure. The frequency of the spectra was estimated ranging from 2 to 150 Hz. Thus, we applied five orthogonal Slepian tapers with TW = 3. The response of spiking activity was defined as the change in mean firing rate relative to the control, by, Rspk =

Rsti − Rblk Rctr

(2)

where Rsti is the mean firing rate of the multi-units under the stimulus, and Rblk is that of the blank stimuli within the same length of time window. The response of LFP was defined as the change in the power spectrum relative to the control (Henrie and Shapley, 2005; Zhang and Li, 2013) by, Plfp =

Psti − Pblk Pblk

(3)

where Psti is the power spectrum (30–65 Hz) of LFP under the stimulus, and Pblk is the power spectrum (30–65 Hz) of the blank stimuli with uniform background, also within the same length of time window. 2.4.2. Size tuning curve The size tuning curve was fitted using the difference of Gaussian (DoG) function (Osaki et al., 2011): R (x) = R0 + Ke × e−(x−e )

2

/2e2

−(x−i )2 /2 2

− Ki × e

i

(4)

where R0 is the spontaneous activity for a blank presentation. Ke ,  e , and e are the amplitude, width, and center of the excitatory component, respectively. And Ki ,  i , and i are the respective equivalents for the surround Gaussian function. Values for Ke ,  e , e , Ki ,  i , and i were optimized to provide least squared error fits to the data. The goodness of fit was assessed using Adjust-R2.

2.4. Data analysis: surround suppression analysis 2.5. Data analysis: estimating connectivity 2.4.1. Individual response Both the firing rate of multi-units and spectra of gamma-band (30–65 Hz) LFP were calculated to determine the surround suppression of individual recording sites. Before the actual analysis, the raw LFP data were preprocessed to minimize the different types of noise (Zhang and Li, 2013). First, the raw LFP for each trial was segmented into small time intervals (10 ms). The segment was considered to be affected by the noise caused by muscles or electrical bursts if it was beyond the mean ± 4 SD, and was then cut. Linear interpolation was subsequently performed to fill in any missing data. The entire trial data set was dropped when more than three noise segments were identified. The 50 Hz line-noise was removed by fitting the LFP to a sine function and subtracting the fitted noise from the signal [Chronux package 2.0 software, (Bokil et al., 2010)]. The spectra of LFP were estimated using the multi-taper method (Percival and Walden, 1993). We used the Chronux package 2.0 software (Bokil et al., 2010), freely download from http://

The functional connectivity was estimated with GC analyses (See Section 2.5) mainly based on gamma-band LFP, because firing rates of some units dropped much lower, underlying the surround modulation conditions, in which case the GC analyses would fail to estimate their causal connection, causing bias for the following analyses. The gamma-band LFP was reconstructed using wavelet transformation. The full-band LFP was first decomposed at five levels with the “db4” wavelet function and the fifth level detail coefficients were then used to reconstruct gamma-band LFPs (31.25–62.5 Hz) also with the “db4” wavelet function. Detailed descriptions and practical operations of GC have been widely reported in previous studies (Barnett and Seth, 2011; Ouyang et al., 2014). In our experiments, we applied the MVGC toolbox to complete the computations. The toolbox could be downloaded freely from www.sussex.ac.uk/sackler provided by Barnett and Seth (2014). The Granger causality between each pair of recording sets

48

X. Niu et al. / Brain Research Bulletin 117 (2015) 45–53

was estimated based on gamma-band LFP pooled data across different trials. The order of the model was 19 and remained the same for all fits in our experiments. After functional connectivity between individual recording sites was estimated, a causal network (indicated with a binary graph, denoted by A) was then transformed by eliminating those interactions not reaching the statistical significance (denoted by ˛) to represent the connection structure of a neuron ensemble. The statistical significance level ˛, namely the threshold, was not kept as a constant in this study, to determine if the conclusion would hold across different levels of significance. We made multiple comparison corrections (False Discovery Rate) when transforming GC matrices to binary graphs across all significance levels. For each estimated causal network, the nodes were defined as the individual recording site, and the edges of the causal network were elements of the binary connection graph. Because the estimated functional connectivity was directed, each node pair of the network may have had up to two edges, one from i to j, and the other one from j to i. Additionally, there were no self-loops, indicating the diagonal elements of the connection matrices were all zero. 2.6. Data analyses: characterizing the network property The graph of the inferred network consisted of N vertices and M directed edges. The averaged connection strength of the population (mGC) and several network parameters (Bullmore and Sporns, 2009) were measured to characterize network properties, including connection density (D), clustering coefficient (C), and average path length (L). 2.6.1. Averaged connection strength Granger causality is asymptotically equivalent to information theoretical transfer entropy (Barnett and Bossomaier, 2012), and often is described as a measure of information flow (Lizier and Prokopenko, 2010). Thus, we calculated mGC by averaging the under-threshold (p < ␣) GC matrix to quantify the amount of information transfer among the networks. The mGC was computed with, 1  gij × wij M N

N

mGC =

(5)

j=1 i=1

where gij indicates the GC value between recording sites I, and j. wij is an element of the binary graph. 2.6.2. Connection density Connection density (represented by D) often indicates the physical cost of a network, for example, the energy or other resource requirements, and is defined as the actual number of edges in the graph as a proportion of the total number of possible edges, which is, M D= 2 N −N

(6)

2.6.3. Clustering coefficient The clustering coefficient (denoted by C) is a measure of the degree to which nodes in a graph tend to cluster together, and was defined as a proportion of the number of connections that existed between the nearest neighbors of a node over the number of all possible connections (Watts and Strogatz, 1998). For the directed network, it was computed with (Fagiolo, 2007),



3

A + AT 1  in+out  in+out ii    C= N 2 K K − 1 − 2 A2 N

i=1

i

i

 ii

(7)

where N is the number of nodes in the network. Kin+out i represents the sum of incoming and outgoing edges at vertex I, and A denotes the binary matrix of the sparse network. 2.6.4. Path length The path length for each pair was defined as the number of edges that formed the shortest path going from one node to the other. The average path length (denoted by L) is inversely related to global efficiency of the network (Achard and Bullmore, 2007), which means the lower value of L, the higher efficiency of parallel information transferred among neurons (Rubinov and Sporns, 2010). The average path length (marked with L) is computed by averaging path lengths across all pairs of nodes (Eq. (8)). L=

 1 di,j 1/2N (N − 1)

(8)

ij

where di , j indicates the shortest distance of the node i and j. 2.7. Histology At the end of each recording, each rat was given an overdose of urethane and perfused transcardially with saline followed by 4% paraformaldehyde in 0.1 M phosphate buffer (pH 7.4). The brain was removed and postfixed for 4–6 h in the same fixative and subsequently cryo-protected by immersion in 30% sucrose in 0.1 M phosphate buffer (pH 7.4) at 4 ◦ C overnight. Each brain was cut coronally with a Leica freezing microtome (CM1950; Leica, Buffalo Grove, IL, USA) at −20 ◦ C to produce 50 ␮m thick sections, which were used to verify the tracks of electrodes (Fig. 1c) using an Olympus inverted fluorescence microscope (IX71; Olympus, Tokyo, Japan). The histological photographs were collected with the Olympus inverted microscope (IX71). Brightness and contrast were adjusted using Adobe Photoshop CS2. 3. Results In our experiments, we completely recorded five neuronal groups from five rat V1 (layer 2/3) areas. The spatiotemporal RF of each neuron was first mapped with sparse noise stimuli (mentioned in Section 2.2) to provide a basis to constitute neuronal groups and further determine the center of the CRF population. The spatial profiles of the RFs at each temporal delay were fitted with 2-D Gabor functions (Eq. (1)). The spatial profile corresponding to the peak response was determined as the final RF profile of the neuron. Because we only took into account individual neurons with clear RF profiles, 14, 13, 14, 13, 15 neurons remained in rats 1–5, respectively. The spatial distributions of RFs from each rat are shown in Fig. 2. The results were analyzed for those five neuron groups. 3.1. Surround modulation of averaged population activity We analyzed the surround modulation of population activity using different sizes of gratings located at RF centers of individual neurons (red dot in Fig. 2). Individual responses were quantified based on spiking activity and LFP. For the spiking activity, the mean firing rate of the multi-units recorded from each electrode was calculated, and the corresponding individual spiking response to each stimulus was computed with Eq. (2) (Section 2.4.1). For LFP signals at the same recording sites, the gamma-band power was calculated and the corresponding individual LFP response was calculated with Eq. (3) (Section 2.4.1). Both responses were then averaged across the population and fitted with the difference of two Gaussian functions (Eq. (4), in Section 2.4.2) to obtain the size-tuning curve for the population.

14

Population averaged LFP gamma-power(μv2)

Population averaged mean firing rate (spks/sec)

X. Niu et al. / Brain Research Bulletin 117 (2015) 45–53

12

10

8

10

20 30 40 Stimulus Diameter (Degree)

a

49

2200

2000

1800

10

20 30 40 Stimulus Diameter(Degree)

b

Fig. 3. Size-tuning of population averaged responses of a sample dataset.

Fig. 3 presents the population-averaged surround modulation effects on both responses based on spiking activity and gammaband LFP for one sample dataset (14 neurons, Rat #1 in Fig. 2). Following the definitions of CRF and eCRF for individual neurons, the grating evoking the strongest population response was treated as CRF (denoted by SCRF ) for the population, and the smallest gratings suppressing response was approximately equal to eCRF (denoted by SeCRF ) for the population. As shown in Fig 3a and b, the surround suppression can be observed clearly in both responses of spiking activity and in the LFP, and the two size tuning curves produced comparable SeCRF (∼40◦ ), and SCRF (∼21.5◦ ). The individual responses modulated as the size of the same set of gratings are shown in Fig. S1a and b, from which we can see that the diminished population activity did not consistently result from feedback from the extra areas to individual neurons. The averaged responses of individual responses described the overall level of the population response. We further did a statistical comparison of SeCRF based on the two kinds of signals for 69 individual neurons, and further compared the results between SCRF based on the two kinds of signals (see Fig. S2). The comparison results indicated that the spiking activity and LFP based results were significantly correlated (r > 0.8, p < 0.001, Spearman’s rank correlation coefficient), for both SCRF and SeCRF . When we calculated the causal connection based on gamma-band LFP (see details in Section 3.2), we mainly relied on the results based on LFP in the subsequent analysis. Thus, for the sample dataset, the surround modulation condition mentioned in the following text is the neuron population with stimulated gratings of a diameter of 40◦ , denoted by CRF+eCRF, while the stimuli denoted by CRF alone indicated gratings of 21.5◦ diameter. The following analyses compared the two stimuli conditions.

3.2. Functional connectivity among neuronal populations The functional connectivity among neuronal populations was estimated with multivariate GC analyses performed upon gammaband LFPs, which were reconstructed with wavelet transforms. The reconstructed gamma-band LFPs of the sample dataset are shown in Fig. 4, from which we can see that the amplitude (power) of the signal under stimulation was much larger than under blank conditions. The estimated GC matrices of the sample datasets for two types of stimuli (CRF alone and CRF+eCRF) are shown in Fig. 5a and e, respectively. The corresponding p-matrices are shown in Fig. 5b and f. A causal network was obtained by discarding those connections not reaching the statistical significance (␣) corrected with

0.5s

Stimulus on Fig. 4. Reconstructed gamma-band LFPs of the sample dataset. The black thick line indicates the time interval under stimulation conditions.

the “False Discovery Rate”. It also serves as a threshold to produce sparse networks, which were more aligned with biological networks. The commonly used significance levels were 0.05 (5%) and 0.01 (1%). In our experiments, the significance level was not kept as a constant to determine whether the variation caused by surround modulation was consistent when the significance level was set at a stricter level. We selected the threshold vector by reducing the ␣ value from 0.01 exponentially to the minimum of a doubleprecision number (˛ = 0.01, 0.001, 0.0001,. . ., 10−16 ). Fig. 5c and g show the connection matrices for the two types of stimuli at ˛ = 0.01, and Fig. 5d and h show connection matrices at ˛ = 10−5 . It was noted that when the ˛ value became much stricter, there was possible emergence of the isolated vertex, not connecting with others in the network. In such a case, we eliminated those and analyzed the properties of the remaining networks. 3.3. Dynamics of functional connectivity modulated by eCRF We analyzed dynamics of functional connectivity modulated by eCRF in this section. Four characteristics of functional connectivity, namely averaged connection strength (mGC), connection density (D), clustering coefficient (C), and path length (L) referred to in Section 2.6, were measured for both stimuli conditions (CRF alone and CRF+eCRF) across different significance levels. Fig. 6 shows

50

X. Niu et al. / Brain Research Bulletin 117 (2015) 45–53

Fig. 5. Transformation from GC matrices to binary connection graphs for two different stimuli conditions of the same orientation.

the statistical comparison results across 12 different orientations between two different stimuli conditions at 16 different significance levels for five rats. To quantify the variation between each group of data, we performed one-way ANOVA analyses at each significance level for each rat. We denoted those not showing significant (0.01 < p < 0.07) variations with arrows (shown in Fig. 6). The untagged groups all showed significant variations (p < 0.01). From Fig. 6 we can see that the variation between two different stimuli conditions for all characteristics were significant even when significance level was set much lower, suggesting that the surround modulation consistently varied the functional connectivity among population networks in V1. As mentioned in Section 2.6, the connection strength (denoted by mGC) was always related to the amount of information transfer among neuronal populations. That is, the large value of mGC indicated more information was transferred among the neuronal population. Furthermore, the connection density (denoted by D) indicated the energy of the networks and the average path length (denoted by L) was inversely related to efficiency of information transfer. The high value of D indicated more causal links in networks, suggesting more energy carried by the network. The low value of L indicated more efficiency of the network to transmit or process parallel visual information. Therefore, the results shown in Fig. 6 indicate that the surround modulation robustly enhanced amount and efficiency of information transfer among the networks across all significance levels.

In addition, the large value of C indicated greater degree to which the neurons tended to cluster together, suggesting the networks were more robust to perturbations (Bullmore and Sporns, 2009). From the comparison results shown in Fig. 6, we can see that cluster coefficients of each network became larger, underlying surround modulation conditions, and suggesting that the network structure was more robust to larger sizes of visual stimuli, which would be beneficial to transfer information for neuron populations.

4. Discussion and conclusions Current studies on surround modulation have emphasized the modulation of individual responses or centers surrounding interactions at the single cell level. In general, these studies have demonstrated that the surrounding regions can modify RF responses through inhibitory effects (Blakemore and Tobin, 1972), spatial summation of low contrast stimuli (Kapadia et al., 1999), cross-orientation modulation (Cavanaugh et al., 2002b; Kimura and Ohzawa, 2008), and even by sharpening the orientation tuning (Okamoto et al., 2009; Liu et al., 2013). In our study, we identified noticeable and robust variations of functional connectivity in V1, resulting from surround modulation upon neuronal populations. The findings generally complemented reports of the neural phenomenon underlying surround modulation at the population level, and provided new insight into the underlying mechanisms.

X. Niu et al. / Brain Research Bulletin 117 (2015) 45–53

51

-3

x 10

11

0.7

0.8

0.6

0.7

0.5

0.6

2.6

10

2.2 2

6

L

D

7

C

mGC

8

0.4

0.5

0.3

0.4

1.8 1.6

5

1.4

4 3

CRF+eCRF CRF alone

2.4

9

0

5

10 |log10(α)|

15

0.2

20

0

5

10 |log10(α)|

15

20

0

5

10 |log10(α)|

15

20

0

5

10 |log10(α)|

15

20

Rat 1 0.02

0.9

0.018

0.95

1.8

0.9

1.7

0.85

1.6

0.8

1.5

CRF+eCRF CRF alone

0.8 0.016 0.7

0.012

L

C

D

mGC

0.014

0.6

0.01

0.75

1.4

0.7

1.3

0.65

1.2

0.5 0.008 0.006

0

5

10 |log10(α)|

15

0.4

20

0

5

10 |log10(α)|

15

20

0

5

10 |log10(α)|

15

1.1

20

0

5

10 |log10(α)|

15

20

Rat 2 -3

11

x 10

0.7

0.8

0.6

0.7

0.5

0.6

2.6

10

2.2

6

2

L

D

7

C

mGC

8

0.4

0.5

0.3

0.4

1.8 1.6

5

1.4

4 3

CRF+eCRF CRF alone

2.4

9

0

5

10 |log10(α)|

15

0.2

20

0

5

10 |log10(α)|

15

20

0

5

10 |log10(α)|

15

20

0

5

10 |log10(α)|

15

20

Rat 3 -3

x 10

0.9

2.6 CRF+eCRF CRF alone

0.8

11

0.8

2.4

10

0.7

0.7

2.2

0.6

2

L

8

C

0.6

D

mGC

9

0.5

0.5

1.8

6

0.4

0.4

1.6

5

0.3

4

0.2

7

0

5

10 |log10(α)|

15

1.4

0.3 0

5

10 |log10(α)|

15

20

0

5

10 |log10(α)|

15

0

5

10 |log10(α)|

15

20

Rat 4 -3

8

2.6

0.8

2.4

0.7

0.7

2.2

0.6

0.6

D

mGC

10

0.9

0.8

6

2

0.5

0.5

1.8

0.4

0.4

1.6

4

0.3

0.3

1.4

2

0.2

0.2

0

5

10 |log10(α)|

15

20

0

5

10 |log10(α)|

15

20

CRF+eCRF CRF alone

L

12

0.9

C

14

x 10

0

5

10 |log10(α)|

15

20

0

5

10 |log10(α)|

15

20

Rat 5 Fig. 6. The variations of characteristics of functional connectivity under two stimuli conditions at different significance levels for 5 rats. Each figure represents the comparison results between the two stimuli conditions, with CRF+eCRF represented with red lines and CRF alone represented with black lines. mGC refers to averaged connection strengths, D connection density, C clustering coefficient, and L path length mentioned in Section 2.6.

52

X. Niu et al. / Brain Research Bulletin 117 (2015) 45–53

Other studies focused on the internal mechanism for organizing or determining surround modulation, and concluded that horizontal connections intrinsic to V1 were likely a key component of the circuit mediating surround suppression (Stettler et al., 2002). In our study, we found that the connection structure among neuron populations in layer 2/3 of V1 where frequently horizontal connections (Burkhalter, 1989; Stettler et al., 2002; Medini, 2011), and varied after surround modulation. These results were consistent with previous reports. In addition, there were studies on the population mapping of contextual modulation, implying that the V1 must be operating as an inhibition-stabilized network (Tsodyks et al., 1997; Okamoto et al., 2009). From our experimental results, we can hypothesize that the V1 may be organized as a dynamic neuronal network, and the spatial pattern of the neuronal network may adapt to different visual tasks. The dynamic functional connectivity has also been found during different phases of working memory tasks (Ouyang et al., 2014), and under different brain states (Nicolaou et al., 2012). Different network activities were also found during shifts in learning strategies (Arias et al., 2014). A static network hypothesis is generally not sufficient to represent the diverse communication patterns of the brain and is therefore not reliable. To estimate the functional or effective connectivity between individual neurons and between areas, one of the most commonly used and powerful techniques is the Granger causality, which could be applied on continuous spiking activity (Ouyang et al., 2014), LFP (Nakhnikian et al., 2014), or even EEG signals (Nicolaou et al., 2012). The method was often combined with graph theory to reveal the characteristics of functional connectivity (Ouyang et al., 2014). Motivated by those previous studies, we quantified the functional connectivity among neuron populations in a similar way. However, we performed further analyses by varying significance levels to determine whether the conclusion would hold when significance level was kept much stricter. The comparison results indicated that the surround modulation did consistently vary the network connectivity activity across all significance levels. Furthermore, the surround modulation conditions for the population were determined according to the averaged value of individual response. The suppression phenomenon of individual response may due to top-down signaling (Wang et al., 2010), although they were not affected by feedback from eCRF (Fig. S1), consistently. The top–down modulation of lateral interactions was also reported to exist in visual cortex (Ramalingam et al., 2013). Therefore, one possible mechanism could explain changes in functional connectivity among V1 neurons, caused by top-down influences. In conclusion, we found that the surround modulation robustly changed the functional connectivity among the neuronal populations in V1, regardless of the strictness of the significance level derived from statistical methods. This study provided additional evidence to understand the neuronal mechanisms of surround modulation at the population level. In general, we hypothesize that the V1 may be organizing as a dynamic neuronal network during the modulation of surround stimuli. Under larger diameter images, neuron groups enhance their efficiency to process or transfer information by increasing their connection density and reducing the path length of the network. However, the conclusions were all obtained from anesthetized animals. It was suggested that the surround or contextual modulation may be mitigated under anesthesia (Vaiceliunaite et al., 2013). However, urethane, a popular anesthetic for acute neurophysiology in rodents (e.g., Sun and Dan, 2009; Zhu and Yao, 2012) did not change synaptic transmission or affect tuning properties (Sceniak and Maciver, 2006). Therefore, anesthesia might not introduce much bias to our analysis. Nevertheless, in future work we will attempt to investigate this issue by using multi-electrode recordings from non-anesthetized animals, which would be more

meaningful for understanding population encoding mechanisms underlying surround modulation or other stimuli conditions. Acknowledgments This study was supported by a Grant (U1304602) from The National Natural Science Funds of China; a Grant (122102210102) from Program for Science and Technology of Henan Province of China; Plan for Scientific Innovation Talents of Henan Province under Grand no. 124200510016. We also thank Dr. Xueguo Zhang (Wayne State University School of Medicine) for his technical assistance in taking the histological photographs shown in Fig. 1c. Appendix A. Supplementary data Supplementary data associated with this article can be found, in the online version, at http://dx.doi.org/10.1016/j.brainresbull. 2015.07.003 References Achard, S., Salvador, R., Whitcher, B., Suckling, J., Bullmore, E., 2006. A resilient, low-frequency, small-world human brain functional network with highly connected association cortical hubs. J. Neurosci. 26, 63–72. Achard, S., Bullmore, E., 2007. Efficiency and cost of economical brain functional networks. PLoS Comput. Biol. 3, e17. Albright, T.D., Stoner, G.R., 2002. Contextual influences on visual processing. Annu. Rev. Neurosci. 25, 339–379. Angelucci, A., Levitt, J.B., Walton, E.J., Hupe, J.M., Bullier, J., Lund, J.S., 2002. Circuits for local and global signal integration in primary visual cortex. J. Neurosci. 22, 8633–8646. Arias, N., Fidalgo, C., Vallejo, G., Arias, J.L., 2014. Brain network function during shifts in learning strategies in portal hypertension animals. Brain Res. Bull. 104, 52–59. Bair, W., Cavanaugh, J.R., Movshon, J.A., 2003. Time course and time-distance relationships for surround suppression in macaque V1 neurons. J. Neurosci. 23, 7690–7701. Barlow, H.B., 1953. Summation and inhibition in the frog’s retina. J. Physiol. 119, 69–88. Barnett, L., Seth, A.K., 2011. Behaviour of Granger causality under filtering: theoretical invariance and practical application. J. Neurosci. Methods 201, 404–419. Barnett, L., Bossomaier, T., 2012. Transfer entropy as a log-likelihood ratio. Phys. Rev. Lett. 109, 138105. Barnett, L., Seth, A.K., 2014. The MVGC multivariate Granger causality toolbox: a new approach to Granger-causal inference. J. Neurosci. Methods 223, 50–68. Blakemore, C., Tobin, E.A., 1972. Lateral inhibition between orientation detectors in the cat’s visual cortex. Exp. Brain Res. 15, 439–40. Bokil, H., Andrews, P., Kulkarni, J.E., Mehta, S., Mitra, P.P., 2010. Chronux: a platform for analyzing neural signals. J. Neurosci.Methods 192, 146–151. Brown, S.P., Hestrin, S., 2009. Intracortical circuits of pyramidal neurons reflect their long-range axonal targets. Nature 457, 1133–1136. Bullmore, E., Sporns, O., 2009. Complex brain networks: graph theoretical analysis of structural and functional systems. Nat. Rev. Neurosci. 10, 186–198. Burkhalter, A., 1989. Intrinsic connections of rat primary visual cortex: laminar organization of axonal projections. J. Comp. Neurol. 279, 171–186. Cavanaugh, J.R., Bair, W., Movshon, J.A., 2002a. Selectivity and spatial distribution of signals from the receptive field surround in macaque V1 neurons. J. Neurophysiol. 88, 2547–2556. Cavanaugh, J.R., Bair, W., Movshon, J.A., 2002b. Nature and interaction of signals from the receptive field center and surround in macaque V1 neurons. J. Neurophysiol. 88, 2530–2546. Chen, G., Dan, Y., Li, C.Y., 2005. Stimulation of non-classical receptive field enhances orientation selectivity in the cat. J. Physiol. 564, 233–243. DeAngelis, G.C., Freeman, R.D., Ohzawa, I., 1994. Length and width tuning of neurons in the cat’s primary visual cortex. J. Neurophysiol. 71, 347–374. Dragoi, V., Sur, M., 2000. Dynamic properties of recurrent inhibition in primary visual cortex: contrast and orientation dependence of contextual effects. J. Neurophysiol. 83, 1019–1030. Fagiolo, G., 2007. Clustering in complex directed networks. Phys. Rev. E Stat. Nonlin. Soft Matter Phys. 76, 26107. Ge, T., Kendrick, K.M., Feng, J., 2009. A novel extended Granger causal model approach demonstrates brain hemispheric differences during face recognition learning. PLoS Comput. Biol. 5, e1000570. Gieselmann, M.A., Thiele, A., 2008. Comparison of spatial integration and surround suppression characteristics in spiking activity and the local field potential in macaque V1. Eur. J. Neurosci. 28, 447–459. Guo, S., Wu, J., Ding, M., Feng, J., 2008. Uncovering interactions in the frequency domain. PLoS Comput. Biol. 4, e1000087.

X. Niu et al. / Brain Research Bulletin 117 (2015) 45–53 Henrie, J.A., Shapley, R., 2005. LFP power spectra in V1 cortex: the graded effect of stimulus contrast. J. Neurophysiol. 94, 479–490. Henry, C.A., Joshi, S., Xing, D., Shapley, R.M., Hawken, M.J., 2013. Functional characterization of the extraclassical receptive field in macaque V1: contrast, orientation, and temporal dynamics. J. Neurosci. 33, 6230–6242. Hubel, D.H., Wiesel, T.N., 1965. Receptive fields and functional architecture in two nonstriate visual areas (18 and 19) of the cat. J. Neurophysiol. 28, 229–289. Hupé, J.M., James, A.C., Girard, P., Bullier, J., 2001. Response modulations by static texture surround in area V1 of the macaque monkey do not depend on feedback connections from V2. J. Neurophysiol. 85, 146–163. Ide, J.S., Zhang, S., Li, C.S., 2014. Bayesian network models in brain functional connectivity analysis. Int. J. Approx. Reason. 56. Jones, J.P., Palmer, L.A., 1987. The two-dimensional spatial structure of simple receptive fields in cat striate cortex. J. Neurophysiol. 58, 1187–1211. Kapadia, M.K., Westheimer, G., Gilbert, C.D., 1999. Dynamics of spatial summation in primary visual cortex of alert monkeys. Proc. Natl. Acad. Sci. U. S. A. 96, 12073–12078. Kimura, R., Ohzawa, I., 2008. Time course of cross-orientation suppression in the early visual cortex. J. Neurophysiol. 101, 1463–1479. Lee, L., Friston, K., Horwitz, B., 2006. Large-scale neural models and dynamic causal modelling. Neuroimage 30, 1243–1254. Liu, K., Yao, H., 2014. Contrast-dependent OFF-dominance in cat primary visual cortex facilitates discrimination of stimuli with natural contrast statistics. Eur. J. Neurosci. 39, 2060–2070. Liu, Y.J., Hashemi-Nezhad, M., Lyon, D.C., 2013. Sharper orientation tuning of the extraclassical suppressive-surround due to a neuron’s location in the V1 orientation map emerges late in time. Neuroscience 229, 100–117. Lizier, J.T., Prokopenko, M., 2010. Differentiating information transfer and causal effect. Eur. Phys. J. B 73, 605–615. Medini, P., 2011. Cell-type-specific sub- and suprathreshold receptive fields of layer 4 and layer 2/3 pyramids in rat primary visual cortex. Neuroscience 190, 112–126. Nakhnikian, A., Rebec, G.V., Grasse, L.M., Dwiel, L.L., Shimono, M., Beggs, J.M., 2014. Behavior modulates effective connectivity between cortex and striatum. PLoS One 9, e89443. Nassi, J.J., Lomber, S.G., Born, R.T., 2013. Corticocortical feedback contributes to surround suppression in V1 of the alert primate. J. Neurosci. 33, 8504–8517. Neapolitan, R.E., 2009. Probabilistic Methods for Bioinformatics: with an Introduction to Bayesian Networks. Foundations of Bayesian Networks. Morgan Kaufmann Publishers, Burlington, MA, pp. 85–134. Nicolaou, N., Hourris, S., Alexandrou, P., Georgiou, J., 2012. EEG-based automatic classification of ‘awake’ versus ‘anesthetized’ state in general anesthesia using Granger causality. PLoS One 7, e33869. Okamoto, M., Naito, T., Sadakane, O., Osaki, H., Sato, H., 2009. Surround suppression sharpens orientation tuning in the cat primary visual cortex. Eur. J. Neurosci. 29, 1035–1046. Okatan, M., Wilson, M.A., Brown, E.N., 2005. Analyzing functional connectivity using a network likelihood model of ensemble neural spiking activity. Neural. Comput. 17, 1927–1961. Osaki, H., Naito, T., Sadakane, O., Okamoto, M., Sato, H., 2011. Surround suppression by high spatial frequency stimuli in the cat primary visual cortex. Eur. J. Neurosci. 33, 923–932. Ouyang, M., Li, S., Tian, X., 2014. Functional connectivity among spikes in low dimensional space during working memory task in rat. PLoS One 9, e91481. Percival, D.B., Walden, A.T., 1993. Spectral Analysis for Physical Applications: Multitaper and Conventional Univariate Techniques. Cambridge Univ., Cambridge, UK, pp. 331–375.

53

Ramalingam, N., McManus, J., Li, W., Gilbert, C., 2013. Top-down modulation of lateral interactions in visual cortex. J. Neurosci. 33, 1773–1789. Reynaud, A., Masson, G.S., Chavane, F., 2012. Dynamics of local input normalization result from balanced short- and long-range intracortical interactions in area V1. J. Neurosci. 32, 12558–12569. Rubinov, M., Sporns, O., 2010. Complex network measures of brain connectivity: uses and interpretations. Neuroimage 52, 1059–1069. Sceniak, M.P., Hawken, M.J., Shapley, R., 2002. Contrast-dependent changes in spatial frequency tuning of macaque V1 neurons: effects of a changing receptive field size. J. Neurophysiol. 88, 1363–1373. Sceniak, M.P., Maciver, M.B., 2006. Cellular actions of urethane on rat visual cortical neurons in vitro. J. Neurophysiol. 95, 3865–3874. Seriès, P., Lorenceau, J., Frégnac, Y., 2003. The silent surround of V1 receptive fields: theory and experiments. J. Physiol. Paris 97, 453–474. Stettler, D.D., Das, A., Bennett, J., Gilbert, C.D., 2002. Lateral connectivity and contextual interactions in macaque primary visual cortex. Neuron 36, 739–750. Sun, W., Dan, Y., 2009. Layer-specific network oscillation and spatiotemporal receptive field in the visual cortex. Proc. Natl. Acad. Sci. U. S. A. 106, 17986–17991. Tsodyks, M.V., Skaggs, W.E., Sejnowski, T.J., McNaughton, B.L., 1997. Paradoxical effects of external modulation of inhibitory interneurons. J. Neurosci. 17, 4382–4388. Vaiceliunaite, A., Erisken, S., Franzen, F., Katzner, S., Busse, L., 2013. Spatial integration in mouse primary visual cortex. J. Neurophysiol. 110, 964–972. Wang, C., Huang, J.Y., Bardy, C., FitzGibbon, T., Dreher, B., 2010. Influence of ‘feedback’ signals on spatial integration in receptive fields of cat area 17 neurons. Brain Res. 1328, 34–48. Watts, D.J., Strogatz, S.H., 1998. Collective dynamics of ‘small-world’ networks. Nature 393, 440–442. Webb, B.S., Dhruv, N.T., Solomon, S.G., Tailby, C., Lennie, P., 2005. Early and late mechanisms of surround suppression in striate cortex of macaque. J. Neurosci. 25, 11666–11675. Wielaard, J., Sajda, P., 2005. Extraclassical receptive field phenomena and short-range connectivity in V1. Cereb. Cortex 16, 1531–1545. Xing, D., Shapley, R.M., Hawken, M.J., Ringach, D.L., 2005. Effect of stimulus size on the dynamics of orientation selectivity in Macaque V1. J. Neurophysiol. 94, 799–812. Yook, C., Druckmann, S., Kim, J., 2013. Mapping mammalian synaptic connectivity. Cell Mol. Life Sci. 70, 4747–4757. Yu, S., Huang, D., Singer, W., Nikolic, D., 2008. A small world of neuronal synchrony. Cereb. Cortex 18, 2891–2901. Zhang, J., Zhang, X., 2010. Electrical stimulation of the dorsal cochlear nucleus induces hearing in rats. Brain Res. 1311, 37–50. Zhang, L., Li, B., 2013. Surround modulation characteristics of local field potential and spiking activity in primary visual cortex of cat. PLoS One 8, e64492. Zhou, D., Xiao, Y., Zhang, Y., Xu, Z., Cai, D., 2014. Granger causality network reconstruction of conductance-based integrate-and-fire neuronal systems. PLoS One 9, e87636. Zhu, Y., Yao, H., 2012. Modification of visual cortical receptive field induced by natural stimuli. Cereb. Cortex 23, 1923–1932.