# Functional connectivity changes during a working memory task in rat via NMF analysis

^{1}School of Biomedical Engineering, Tianjin Medical University, Tianjin, China^{2}Research Center of Basic Medicine, Tianjin Medical University, Tianjin, China

Working memory (WM) is necessary in higher cognition. The brain as a complex network is formed by interconnections among neurons. Connectivity results in neural dynamics to support cognition. The first aim is to investigate connectivity dynamics in medial prefrontal cortex (mPFC) networks during WM. As brain neural activity is sparse, the second aim is to find the intrinsic connectivity property in a feature space. Using multi-channel electrode recording techniques, spikes were simultaneously obtained from mPFC of rats that performed a Y-maze WM task. Continuous time series converted from spikes were embedded in a low-dimensional space by non-negative matrix factorization (NMF). mPFC network in original space was constructed by measuring connections among neurons. And the same network in NMF space was constructed by computing connectivity values between the extracted NMF components. Causal density (Cd) and global efficiency (E) were estimated to present the network property. The results showed that Cd and E significantly peaked in the interval right before the maze choice point in correct trials. However, the increase did not emerge in error trials. Additionally, Cd and E in two spaces displayed similar trends in correct trials. The difference was that the measures in NMF space were significantly greater than those in original space. Our findings indicated that the anticipatory changes in mPFC networks may have an effect on future WM behavioral choices. Moreover, the NMF analysis achieves a better characterization for a brain network.

## Introduction

Working memory (WM) has been defined as a limited capacity system which temporarily maintains and manipulates information “online” during various cognitive tasks (Baddeley, 2003). It supports human thought processes by providing an interface between perception, long-term memory and action. Consequently, the neural activity supporting WM underlies a broad range of advanced cognitive functions (Jonides et al., 2008).

The Y-maze spontaneous alternation task is one of the most widely used WM paradigms. Baeg et al. (2003) have found that firing rate changes of PFC neurons is tied to the reward and reward expectancy during a spatial delayed alternation task on a figure 8-shaped maze, implicating that WM is mediated by the PFC neuronal activities. Benchenane et al. (2010) have suggested that prefrontal cortex (PFC) is involved in rats learning new rules on a Y-maze both at the local field potential and spike levels. There is also evidence that PFC neuronal activity displays anticipatory changes prior to future goal-directed behaviors, suggesting that the PFC is responsible for holding neural representations for WM (Jo et al., 2013).

During a cognitive process, neuronal discharges are not independent but serially correlated with each other (Constantinidis and Goldman-Rakic, 2002). Therefore, a brain can be functionally defined as a network. Studies of functional network organization have become one of the most popular topics in the field of neuroscience (Palva et al., 2010; Sakaki et al., 2013). Several research has shown that cognitive performance was mediated by interactions among distributed, functionally specialized brain networks (White et al., 2011; Sala-Llonch et al., 2012). Conversely, specific network disruptions can cause a range of neurological disorders (Seeley et al., 2009). A powerful technique called Granger causality analysis could be used to extract such connections. Brain connectivity in humans or animals during cognitions has been analyzed using EEG (Anderson et al., 2010), fMRI (Schlösser et al., 2006), MEG (Bressler and Seth, 2011) and LFPs (Bressler et al., 2007), providing satisfactory explanations for neural mechanisms. These studies have opened new thought that the functional connectivity in mPFC networks may play an essential role in the WM. With the rapid development of chronic multi-electrode implantation, spikes can be recorded extracellularly. It is well-known that the spikes have better spatial and temporal resolutions than other signals. Therefore, these advances in the electrode implantation and brain network analyses have allowed us to map brain functional connectivity patterns *in vivo*.

A previous study has demonstrated that selective interaction among neurons is an essential component in neural information processing (Benchenane et al., 2010). Based on Hebb's theories, only a small fraction of neurons (“cell assembly”) participate in a task, selecting and preserving most valuable information to conduct future behaviors. Previous studies have demonstrated that sparse-coding is a convincing hypothesis to explain the way to process neural information (Sakurai, 1999; Schweighofer et al., 2001). Several studies have also reported that brain networks are sparse (Sporns et al., 2005; Hwang et al., 2013). However, most studies have only focused on functional connectivity among all recorded neurons (Lin et al., 2009; Nedungadi et al., 2009; Kim et al., 2011). It is cautioned that spurious effects may be introduced by neurons not in assemblies. From a network perspective, it is necessary to separate present vs. absent connections. Many investigators have attempted to extract the true network features by several mathematics methods. Recent work shows that cell assemblies can be extracted, which is used to represent the intrinsic pattern of neuronal activity (Peyrache et al., 2009; Benchenane et al., 2010). Up to now, non-negative matrix factorization (NMF), principal component analysis (PCA) and independent component analysis (ICA) have been already used to extract the latent structure of fMRI data (Lohmann et al., 2007; Shen and Meyer, 2008; Demirci et al., 2009; Zhou et al., 2009). In our study, NMF is employed to learn the parts-based components of realistic objects with non-negative constraints (Lee and Seung, 1999). The whole perception in brain relies on its parts-based representations (Logothetis and Sheinberg, 1996; Lee and Seung, 1999). And neuronal firing rates are never negative. Based on the theories, NMF has found its way into the neuroscience community (Lohmann et al., 2007; Padilla et al., 2010; Liu et al., 2013). It can be seen that the notion of combining parts to form a whole is intuitive and easy-to- interpreted. Here, NMF was used to remove dimensions that did not contribute to the raw firing rate pattern of mPFC neurons. The necessary information was correctly retained by NMF components. Therefore, the method could not only find intrinsic feature of neural activity, but also reduce computational cost.

In summary, a deeper understanding of the mPFC network changes may have profound implications for revealing the neural mechanism underlying WM. In this present work, we examined connectivity dynamics in mPFC spike networks in both NMF and original spaces, in order to investigate: (1) the association between the network changes and WM behavioral choices; (2) whether the network feature can be better represented in NMF space.

## Experimental Procedures

### Subjects

Five male Sprague-Dawley rats weighting 300–350 g participated in this study. They were housed in temperature-controlled and clear plastic cages (25°C) under a 12:12-h light-dark schedule. All Experiments were conducted during the light phase. The experimental protocol was approved by the Committee on Experimental Animal of Tianjin Medical University.

### Apparatus

The apparatus was a Y-maze (Figure 1A), consisting of a start arm A and two goal arms B and C. Three arms were identical in size (80 cm long, 16 cm wide and 21 cm high) and other physical features. They were linked together in a symmetrical Y shape (120° of angular deviation from each other). A starting platform on the arm A was separated by a guillotine door. Two food wells were located at the end of each goal arm, respectively. Meanwhile, there was an infrared detector at the fork of the Y-maze to record a reference time (RT). The RT referred to the time for choice behavior.

**Figure 1. Working memory behavior**. **(A)** A working memory task on a Y-maze (overhead view). An entire trial contained a free choice run (left) and an alternative choice run (right). Left: Arm A is a start arm. Reward sites are located at the ends of two goal arm B and C. RT was recorded by an infrared detector (red triangle). Right: A black line represents a correct choice and a gray line represents an error choice. **(B)** As training processed, the average rate of correct choice gradually increased (*n* = 5 rats). **(C)** As training continued, the average latency (*n* = 5 rats) gradually decreased. **(D)** The latency distribution (from the guillotine door to the choice point, *n*_{Rat1} = 24 trials, *n*_{Rat2} = 15 trials, *n*_{Rat3} = 21 trials, *n*_{Rat4} = 32 trials, *n*_{Rat5} = 22 trials) for each rat. **(E)** Average latency (from the guillotine door to the reward site, *n*_{Rat1} = 24 trials, *n*_{Rat2} = 15 trials, *n*_{Rat 3} = 21 trials, *n*_{Rat4} = 32 trials, *n*_{Rat 5} = 22 trials) over entire trials for each rat.

### Behavioral Experiments

The behavioral task was a spatial alternation Y-maze paradigm. First, all the rats were provided water *ad libitum* and food restriction for 2 days, maintaining at 85% of free-feeding body weight. Next, the rats explored the Y-maze and consumed reward food over two consecutive days (30 min per day). During habituation, they were forced to perform at least 15 entries to reward arms in 30 min. After that, the rats were trained to learn an alternative choice rule on the Y-maze. They were given one training session of 10 entire trials per day. An entire trial consisted of a free choice run and an alternative choice run. A latency for each trial referred to a period when a rat started from start point to reward site. During the free run, entering either goal arms was all rewarded. After ate, the rats ran back to the start location spontaneously and a guillotine door was closed. After 5 s delay, the door was lifted and the rats started to perform the choice run. They must choose a goal arm not previously entered to get reward. An error was defined by re-entering the arm visited in the last free run. On the way home, the arm that had not been visited was blocked by a guillotine door in order to prevent the rats quickly turning around and running into this arm. When the percent correct reached at least 80% over consecutive 2 days, the rats were considered as acquiring the rule used for the Y-maze task. If the rat has not learned the rule on day 8, it would be knocked out.

### Electrode Implantation

After reaching the criterion of 80%, the rats were deeply anesthetized with sodium pentobarbital (40 mg/kg) supplemented with atropine (0.2 mg/kg) injected abdominally. Then a microelectrode array (2 × 0.3 mm2 in area, made in house) was chronically implanted in the medial PFC (mPFC, [AP: 2.5–4.5 mm; ML: 0.2–1 mm; DV: 2.5–3.5 mm]). This array consisted of 16 wires (2 × 8, 33 micron diameter nickel-cadmium, spaced 200 micron, <1 MΩ). During surgical procedures, animal's body temperature was monitored and maintained at 35°C~37°C. After the surgery, the rats were allowed to recover for at least 2 weeks. Food and water were available *ad libitum*.

### Electrophysiological Recording and Data Preprocessing

Following 2 weeks recovery, the rats were again conducted into performed the Y-maze task. At the same time, neuronal activity was recorded. The RT for each trial was detected by an infrared sensor mounted on the Y-maze. The latencies for all trials were also recorded. To be clear, the trials mentioned in following paper referred to the alternative choice trials which could be used to measure the WM performance. A total 114 trials were analyzed (*n*_{correct in total} = 96, *n*_{error in total} = 18. Rat 1: *n*_{correct} = 20, *n*_{error} = 4; Rat 2: *n*_{correct} = 13, *n*_{error} = 2; Rat 3: *n*_{correct} = 18, *n*_{error} = 3; Rat 4: *n*_{correct} = 27, *n*_{error} = 5; Rat 5: *n*_{correct} = 18, *n*_{error} = 4) in the present paper. Besides, the baseline activity was taken during the inter-trial intervals. We compared the certain duration of mPFC neuronal activity (duration: 0.25 s pre—measure peak and 0.25 s post—measure peak) in correct trials with the baseline activity.

During the recording courses, the rats were connected to the Cerebus Acquisition system (Cyberkinetics, Foxborough, MA) via a rotating and lightweight cable. They were allowed to move freely on the Y-maze. Wide-band neural signals through the multi-electrode array were fed into a head stage amplifier with 5000 × gain and band-pass filtered between 250 and 7.5 kHz. The obtained signals whose amplitudes exceeded pre-set voltage thresholds were recorded, digitized at 30 kHz and stored for off-line spike sorting. Pre-set voltage thresholds for spike detection were done on negative amplitudes for each electrode independently. The threshold level was set as –7.5 times the standard deviation of raw signals (Reyes-Puerta et al., 2014). Any waveforms that crossed the corresponding threshold were detected and stored with time stamps. Since neural activity per channel was actually synthesized by several neurons nearby an electrode tip, Off-Line Sorter software (Plexon) was applied to extract multiple spike trains associated with different neurons by a semiautomatic clustering selection method (Lewicki, 1998) followed by manual clustering.

### Histology

After electro-physiological recordings, all the rats were sacrificed for histological analyses. The animals anaesthetized with sodium pentobarbital (40 mg/kg) were perfused intracardiacally with saline followed by 10% formalin saline. Then the brains were removed from the calvaria and trimming. Serial, 150 μm thick, coronal sections were carried out. Tracks and lesion sites were observed under a light microscope.

### Construction of a Continuous Time Series From a Spike Train

Due to the discrete nature, spikes were not suitable for the subsequent connectivity analysis. A spike train must be converted into a continuous instantaneous firing rate time series (Zhu et al., 2005). The new series could encode the same temporal firing pattern of neuron as original spikes.

Specifically, the mathematical procedure consisted of three steps (Zhu et al., 2005). First, the instantaneous firing rate could be described by the following function:

where *t*_{i − 1}, *t _{i}* and

*t*

_{i + 1}were considered as three successive spike times. Next, a small time interval δ

*T*was chosen to obtain the instantaneous integrated rate as follow:

Note that δ*T* was much smaller than the mean interspike interval *T*. According to the recommendation (Zhu et al., 2005), δ was set to 0.25. Finally, a continuous instantaneous firing rate series was yielded by smoothing out above functions in time (by Lagrange interpolation methods). In this way, each trial was represented by a matrix which consisted of *m* neurons of *n* instantaneous firing rates.

### Extracted Components in NMF Space

*X* is an *m*-by-*n* input matrix. *m* is the number of neurons. *n* is the number of instantaneous firing rates. NMF was performed on the input matrix to decompose it into non-negative basis matrix (*W*) and non-negative encoding matrix (*H*), as follows:

where FRO is known as Frobenius norm and *e* is residual error matrix (Paatero and Tapper, 1994; Lee and Seung, 1999). *W* is a representation of NMF space. *H* is a low-dimensional representation of the input data (Padilla et al., 2010).

In NMF analysis, the estimation of rank *r* for factorization is a critical matter. It is required that at least:

The NMF computation was iterated by varying the parameter *r* from 1 to *m* (the number of neurons). The number whose variance accounted for (*VAF*) was just greater than 90% was selected as the rank *r* (Torres-Oviedo et al., 2006). The *VAF* is expressed as follows:

For each trial, there was a NMF space. The *VAF* and rank *r* for each one were calculated. The rank *r* varied between animals, since different numbers of neurons were recorded as well as neurons with different contributions. There were also different ranks constructed from the same animal in different trials, as different neurons were either active or silent in different trials, or patterns of activity changed. In order to perform statistical analysis on data, maximum rank *r* for each rat was determined as the number of NMF components which keeps the *VAF* for each trial more than 90%. That means the neural activity pattern was well accounted for by extracted NMF components. Thus, the rank *r* was consistent within an animal, across trials. According to the final rank, we recalculated the corresponding *VAF* for each trial.

### Network Analysis

#### Network construction

Granger causal connectivity analysis (GCCA) is a powerful approach for assessing functional connectivity among neural signals, revealing how one time series might exert influence on another (Zhou et al., 2013). If the past of one series does not contribute to the future of another one, the Granger causality (GC) between them is 0. When applied to a brain, this method has a potential to map the connectivity pattern of a brain network. A network is defined as a set of nodes linked by edges. In original space, nodes represent a set of neurons recorded in mPFC. In NMF space, nodes represent a group of NMF components. Edges represent the Granger causality values between their respective nodes. A 0.5-s moving window (moving step: 0.125 s) was applied to time series of every trials. Here we give a brief mathematical description about the Granger causality analysis (Seth, 2010). First, the multivariate time series can be described as follows:

Next, a signal *X* is represented by a multivariate autoregressive model as follows:

where *p* is the model order specified by the Bayesian information criterion, *A(n)* is the model coefficients, and *E(t)* contains the prediction errors for all time series, and ξ_{m}(*t*) is the prediction error for *x _{m}*. The

*E*(

*t*) covariance matrix of the unrestricted model is defined as follows:

In a restricted model omitting variable *x _{j}*,

*E*(

*t*) is represented as follows:

Thus, the *E _{R}*(

*t*) covariance matrix of the restricted model omitting

*x*is:

_{j}The Granger causality from variable *x _{j}* to variable

*x*, conditioned on all other variables is given by:

_{i}Having constructed a Granger causality matrix, it is important to assess the statistical significance. Only significant connections were retained by *F-*tests (Seth, 2010; Frye et al., 2012). The tests were corrected for multiple comparisons by means of the false discovery rate (Benjamini and Hochberg, 1995). Connections which were not statistically significant were set to zeros.

#### Network properties

In the present study, we focus on the global properties in mPFC networks during the WM task. Both causal density (Cd) and global efficiency (E) can provide overall descriptions of a network.

Cd is a measure of causal interactions involving dynamical integration (neurons coordinate with each other) and differentiation (neurons contribute in different ways) among network nodes. Given a causal network identified by GCCA, Cd is defined as (Seth et al., 2006):

where *Num _{GC}* is the total number of significant connections,

*N*is the number of nodes. The values of

*Cd*are bounded between 0 and 1. A neural network only containing independent neurons would have low

*Cd*. By contrast, a network with high

*Cd*indicates that nodes globally coordinate with each other and are dynamically distinct (Seth, 2010) to accomplish a task at the same time.

E is a measure of the capacity for nodes to propagate information in parallel across a network (Hwang et al., 2013). In general, the integration feature of a network is based on paths, which represent potential routes of information flow between pairs of nodes (Rubinov and Sporns, 2010). The path length is the sum of edge lengths along this path. Every edge length is determined by calculating the inverse of connection strength. The shortest path length between node *i* and *j* is the path with minimum total length (Yan et al., 2011). *E* could be defined mathematically as Bernhardt et al. (2011):

where *L _{ij}* is the shortest path length between node

*i*and

*j*,

*N*is the number of nodes. Previous studies had demonstrated that

*E*was more suitable for sparser networks (Achard and Bullmore, 2007; Estrada and Hatano, 2008; Rubinov and Sporns, 2010). As the mPFC network in our study is sparse, inversing path length simplifies the numerical issues in the computation of

*E*in sparser networks. Moreover,

*E*measures the capacity to transfer parallel information. Since the brain supports large amounts of parallel information processing,

*E*is a superior measure of network integration. A network having a high

*E*suggests the existence of strong associations between nodes. Information is communicated with a high level of efficiency throughout the network.

### Statistical Analysis

In this study, all data was expressed as mean ± standard error (SEM). Statistical analysis was carried out using One- or Two-Way ANOVA for repeated measures of the network features. *T*-test was used to compare the measures in the WM and baseline states. The significance level was set at *P* < 0.05.

### NMF-GCCA Algorithm

In order to be more clearly, the general procedure for researching the mPFC network changes in original and NMF spaces was summarized as follows:

(1) Obtain spike trains in rat mPFC during a WM task and convert each discrete train into a continuous time series.

(2) Project multiple continuous series into a low-dimensional space by NMF, yielding NMF components which can represent the intrinsic pattern of mPFC neurons in original space.

(3) Construct mPFC networks in NMF and original spaces using GCCA.

(4) Analyze the properties of mPFC networks by computing *Cd* and *E*.

(5) Repeat above steps for all trials of 5 rats.

(6) Comparison of mPFC network dynamics during correct trials in both NMF and original spaces.

(7) Compare mPFC network features between the certain WM and baseline states.

(8) Compare mPFC network properties in correct and error trials.

## Results

### Working Memory Behavior

As training processed, the average correct rate (*n* = 5 rats) improved from 63 ± 5 to 91 ± 8% over 8 days (Figure 1B). And the average latency (*n* = 5 rats), from the start point to the choice point, gradually fell from 11.17 ± 0.88 s to 1.94 ± 0.12 s (Figure 1C). Finally, 5 rats successfully learned to alternate between the two goal arms to obtain reward (Figure 1A). By this time, the latency distribution (Rat 1: *n* = 24 trials; Rat 2: *n* = 15 trials; Rat 3: *n* = 21 trials; Rat 4: *n* = 32 trials; Rat 5: *n* = 22 trials) for each rat is shown in Figure 1D. The average latencies over a full trial were 1.6 ± 0.09 s, 1.7 ± 0.15 s, 1.5 ± 0.12 s, 1.7 ± 0.11 s, and 2.2 ± 0.25 s, respectively (Figure 1C).

### Neuronal Activity in mPFC During Working Memory Tasks

Spikes in mPFC were recorded while rats performed WM tasks. Under a light microscope, the electrode tracks and lesion sites were unequivocally identified for 5 rats. The region of an implanted electrode array is provided in Figures 2A–D. Figure 3A displays representative spike activity in a correct trial. 15-channel spike trains were sorted, obtaining 24 single-unit spike trains (Figure 3B). Based on sorting results, the relationship between channel numbers and neuron numbers is shown in Figure 3C. For example, two neurons (neuron 13 and 14) were isolated from the Channel 7 (Figure 3D). Under this method, in total of 117 neurons (Rat 1: *n* = 17; Rat 2: *n* = 24; Rat 3: *n* = 32; Rat 4: *n* = 23; Rat 5: *n* = 21) were collected for analysis. The firing rates of correct trials were averaged using a series of sliding windows (width: 0.5 s, moving step: 0.125 s; Rat 1: *n* = 20 trials; Rat 2: *n* = 13 trials; Rat 3: *n* = 18 trials; Rat 4: *n* = 27 trials; Rat 5: *n* = 18 trials). One-Way ANOVA shows main effect of time on average firing rates (Rat 1—Rat 5, *p* < 0.05), indicating that the values significantly peaked before rats arriving at the maze choice point.

**Figure 2. Location of the electrode array in mPFC**. **(A)** Left, a rat brain model (overhead view). Red point indicates the Bregram. 2 × 8 black dots represent a micro-electrode array. **(B)** A representative photograph of a rat brain. The black arrow illustrates the extent of mPFC lesion after an electrode array was moved through this tissue. **(C)** A coronal section of a rat brain. An electrode array (red arrow) was implanted in PrL (black arrow). **(D)** Photomicrograph showing the tracks (red arrows) of recording electrodes.

**Figure 3. Multi-channel spike sorting**. **(A)** Multi-unit spike trains from 15 channels during a correct trial. Time 0 (red vertical line) indicates the RT. **(B)** Single-unit spike trains sorted from the recordings in **(A)**. **(C)** A corresponding relationship between the electrode number and the neuron number. **(D)** Two neurons were isolated from the Channel 7. **(E)** Average firing rate as a function of time in correct trials for 5 rats (*n*_{Rat1} = 20 trials, *n*_{Rat2} = 13 trials, *n*_{Rat3} = 18 trials, *n*_{Rat 4} = 27 trials, *n*_{Rat5} = 18 trials).

### Embedding Dataset in NMF Space

Figure 4A shows single-unit spikes during a WM trial. The spike trains were then converted into continuous-time series via the instantaneous firing rate method. Two representative signals (Neuron 8 and 16) are shown in Figure 4B. The results suggest that continuous series can describe the temporal firing patterns of neurons. Then the resulting time series were embedded in NMF space. Figure 4C depicts that *VAF* varied with rank *r*. By using the criterion previously described, the data shown in Figure 4A has a rank of 7. According to the selected rank, 7 NMF components were identified from this trial (Figure 4D).

**Figure 4. Data preprocessing and NMF component extraction. (A)** The rastergram shows single-unit spike trains during a trial. **(B)** Two examples of constructions of continuous time series from discrete spike trains (Neuron 8 and 16 from **A**). **(C)** The percentages of Variance Accounted For with various number of rank r. Red dashed line means the threshold for *VAF* (90%). 7 (red point) is the optimized number of rank r in this trial. **(D)** NMF components extracted from the continuous series.

In order to keep a same rank across trials in a rat, the maximum rank *r* was chosen to be a final rank *r* (Figure S1). They were 5, 7, 10, 8, and 7 for Rat #1, #2, #3, #4, and #5. The corresponding mean *VAF* was 95.1 ± 0.84% (Rat 1, *n* = 20 trials), 96.3 ± 0.92% (Rat 2, *n* = 13 trials), 91.8 ± 0.26% (Rat 3, *n* = 18 trials), 95.7 ± 0.47% (Rat 4, *n* = 27 trials) and 92.9 ± 0.48% (Rat 5, *n* = 18 trials), respectively. The results show that based on the maximum rank, neural activity was well accounted for by extracted NMF components.

### Changes of mPFC Networks During Correct Trials in Original and NMF Spaces

mPFC networks in both original and NMF spaces were constructed by GCCA method. The mean connectivity matrices (duration: 0.25 s pre—measure peak and 0.25 s post—measure peak) during correct trials in two spaces are shown in Figure 5A (original space) and Figure 5B (NMF space) for 5 rats (Rat 1: *n* = 20 trials; Rat 2: *n* = 13 trials; Rat 3: *n* = 18 trials; Rat 4: *n* = 27 trials; Rat 5: *n* = 18 trials). From visual inspection, the matrices in NMF space appear to exhibit strong connections. The subsequent statistics for results suggest that mean connection values during the certain WM durations in NMF space (Rat 1, 0.15 ± 0.015; Rat 2, 0.12 ± 0.010; Rat 3, 0.15 ± 0.01; Rat 4, 0.14 ± 0.01; Rat 5, 0.15 ± 0.013) are significantly greater than those in original space (*t*-test, Rat 1, 0.04 ± 0.005, *p* < 0.01; Rat 2, 0.04 ± 0.005, *p* < 0.01; Rat 3, 0.04 ± 0.002, *p* < 0.01; Rat 4, 0.04 ± 0.003, *p* < 0.01; Rat 5, 0.06 ± 0.009, *p* < 0.01). Additionally, the mean connection weights as a function of time in correct trials are shown in Figure 5C. It can be found that the measures peak before correct behavioral choices. And the changes of mean connection weights in NMF and original spaces appear to track one another. However, compared with the values in original space, the NMF analysis did increase the average connection weights (Two-Way ANOVA, factors: NMF / original space (S) and time on maze (T); Rat 1, *n* = 20 trials, main effect: S: *p* < 0.05, T: *p* < 0.05, interaction *p* < 0.05; Rat 2, *n* = 13 trials, main effect: S: *p* < 0.01, T: *p* < 0.05, interaction *p* < 0.05; Rat 3, *n* = 18 trials, main effect: S: *p* < 0.01, T: *p* < 0.01, interaction *p* < 0.01; Rat 4, *n* = 27 trials, main effect: S: *p* < 0.01, T: *p* < 0.01, interaction *p* < 0.01; Rat 5, *n* = 18 trials, main effect: S: *p* < 0.01, T: *p* < 0.01, interaction *p* < 0.05), especially near the peaks.

**Figure 5. mPFC spike networks during WM in original and NMF spaces. (A)** Mean connectivity matrices during WM for 5 rats in original space (*n*_{Rat1} = 20 trials, *n*_{Rat 2} = 13 trials, *n*_{Rat 3} = 18 trials, *n*_{Rat 4} = 27 trials, *n*_{Rat5} = 18 trials). Each node represents a neuron. Colored grid represents the Granger causality between each pair of neurons. The hot colors show strong connections and the cold ones show weak connections. **(B)** Mean connection matrices in NMF space (*n*_{Rat 1} = 20 trials, *n*_{Rat 2} = 13 trials, *n*_{Rat 3} = 18 trials, *n*_{Rat 4} = 27 trials, *n*_{Rat 5} = 18 trials). Each node represents a NMF component. Colored grid represents the Granger causality between each pair of NMF components. The grids are coded according to the connection weight. **(C)** Average Granger causality weights as a function of time in correct trials for 5 rats (red: NMF space, blue: original space; *n*_{Rat 1} = 20 trials, *n*_{Rat 2} = 13 trials, *n*_{Rat 3} = 18 trials, *n*_{Rat 4} = 27 trials, *n*_{Rat 5} = 18 trials).

### Causal Density Dynamics of mPFC Networks During Correct Trials in Original and NMF Spaces

Similar to the mean connection weights, Causal density (Cd) increased to peaks in the interval right before rats' correct behavioral choices. The distinct increases were observed in both original and NMF spaces for 5 rats (Figure 6, from left to right; Two-Way ANOVA, factors: NMF / original space (S) and time on maze (T); Rat 1, *n* = 20 trials, main effect: S: *p* < 0.01, T: *p* < 0.05, interaction *p* < 0.05; Rat 2, *n* = 13 trials, main effect: S: *p* < 0.01, T: *p* < 0.01, interaction *p* < 0.05; Rat 3, *n* = 18 trials, main effect: S: *p* < 0.01, T: *p* < 0.01, interaction *p* < 0.01; Rat 4, *n* = 27 trials, main effect: S: *p* < 0.01, T: *p* < 0.01, interaction *p* < 0.01; Rat 5, *n* = 18 trials, main effect: S: *p* < 0.01, T: *p* < 0.05, interaction *p* < 0.05). Additionally, Cd in NMF space was significantly higher than that in original space. The results indicated that the changes of Cd in correct trials were represented more obvious by NMF analysis.

**Figure 6. Causal density (Cd) of mPFC networks during correct WM trials in original and NMF spaces**. Red lines: NMF space. Blue lines: original space (*n*_{Rat 1} = 20 trials, *n*_{Rat 2} = 13 trials, *n*_{Rat 3} = 18 trials, *n*_{Rat 4} = 27 trials, *n*_{Rat 5} = 18 trials). *Cd* values in the two spaces display a similar changing tendency across correct trials. They increased to peaks before the RT (time 0, black vertical line). Cd in NMF space was significantly higher than that in original space. *Abscissa*: time course that the rat ran from the start point to the reward site. Shaded area indicates the SEM.

### Global Efficiency Changes of mPFC Networks During Correct Trials in Original and NMF Spaces

Global efficiency (E) changes of mPFC networks were examined in both original and NMF spaces for each rat (Figure 7). During correct trials, we found the trends of increased E occurring before the choice point in two spaces (Figure 7, from left to right; Two-Way ANOVA, factors: NMF / original space (S) and time on maze (T); Rat 1, *n* = 20 trials, main effect: S: *p* < 0.01, T: *p* < 0.05, interaction *p* < 0.05; Rat 2, *n* = 13 trials, main effect: S: *p* < 0.01, T: *p* < 0.05, interaction *p* < 0.05; Rat 3, *n* = 18 trials, main effect: S: *p* < 0.01, T: *p* < 0.01, interaction *p* < 0.01; Rat 4, *n* = 27 trials, main effect: S: *p* < 0.01, T: *p* < 0.01, interaction *p* < 0.01; Rat 5, *n* = 18 trials, main effect: S: *p* < 0.01, T: *p* < 0.01, interaction *p* < 0.05). The statistical results indicate that there were profound modifications in E of mPFC networks during correct trials. And like the Cd, E in NMF space was significantly greater than that in original space, indicating that the WM-related mPFC network dynamics was highlighted in NMF space.

**Figure 7. Global efficiency (E) features of mPFC spike network during correct WM trials in original and NMF spaces**. Red lines: NMF space. Blue lines: original space (*n*_{Rat 1} = 20 trials, *n*_{Rat 2} = 13 trials, *n*_{Rat 3} = 18 trials, *n*_{Rat 4} = 27 trials, *n*_{Rat 5} = 18 trials). E increased to peaks prior to the maze choice point in both of the two spaces. The difference was that the mean E in NMF space stayed above the values in original space. *Abscissa*: time course that the rat ran from the start point to the reward site. Shaded area indicates the SEM.

### mPFC Network Properties in Certain WM and Baseline States

The above results demonstrated that the Cd and E in mPFC networks significantly increased before the maze choice point during correct trials. In order to examine whether the increases were meaningful changes, we compared the measures in certain periods (duration: 0.25 s pre—measure peak and 0.25 s post—measure peak, *n* = 96 trials, 5 rats) in correct trials to measures in the baseline states (*n* = 91 trials, 5 rats). As shown in Figure 8, *t*-test reveals that there are significant differences between the certain WM and baseline groups (Figure 8A, Cd: *p* < 0.01; Figure 8B, E: *p* < 0.01).

**Figure 8. Causal density (Cd) and global efficiency (E) of mPFC networks in certain WM and baseline states in NMF space**. Red: certain WM states. Black: baseline states. **(A)** Cd of certain WM states in correct trials (duration: 0.25 s pre-measure peak and 0.25 s post-measure peak, *n* = 96 trials, 5 rats) is significantly greater than that in baseline states (*n* = 91 trials, 5 rats, *t*-test, ***p* < 0.01). **(B)** E of certain WM states in correct trials (duration: 0.25 s pre-measure peak and 0.25 s post-measure peak, *n* = 96 trials, 5 rats) is significantly higher than that in baseline states (*n* = 91 trials, 5 rats, *t*-test, ***p* < 0.01).

### mPFC Network Properties in Correct and Error WM Trials

Since Cd and E in mPFC networks changed dynamically during correct trials, we attempted to ascertain the network properties when rats made error choices. The results show that measures in error groups are significantly lower than those in correct groups (Figures 9A,B; *n*_{Correct} = 96 trials; *n*_{Error} = 18 trials, Two-Way ANOVA, factors: correct/error choice (C) and time on maze (T); Cd, main effects: C: *p* < 0.01, T: *p* < 0.01, interaction *p* > 0.05; E, main effects: C: *p* < 0.01, T: *p* < 0.01, interaction *p* > 0.05; 5 rats). And the Cd and E almost stayed constant with slight oscillation over entire time courses. One-Way ANOVA was then applied to determine whether there was a main effect of time on the measures in error trials. The results indicate that no distinct changes were found in error trials (*n*_{Error} = 18 trials, One-Way ANOVA: Cd: *p* > 0.05; E: *p* > 0.05; 5 rats).

**Figure 9. mPFC network properties in correct and error trials in NMF space**. Red: correct trials. Green: error trials. **(A)** Mean Cd in correct (*n* = 96, 5 rats) vs. error (*n* = 18, 5 rats) trials. No significant increase was observed in error trials. *Abscissa*: the minimum average latency among 5 rats. Shaded area indicates the SEM. **(B)** Average E in correct (*n* = 96, 5 rats) vs. error (*n* = 18, 5 rats) trials. The values showed a slight concussion across the whole trials. *Abscissa* and shaded area: Same as in **(A)**.

## Discussion

This present study explored mPFC network changes during a WM task, employing multi-channel spikes. We embedded spike trains into a NMF space, followed by network analysis in this transformed space. The principle findings here are (1) Cd and E of mPFC networks changed dynamically while the rats correctly performed WM tasks; (2) before the rats arrived at the maze choice point, Cd and E had already augmented to peaks; (3) mPFC network properties were highlighted by NMF analysis. Furthermore, mPFC network measures of certain WM states in correct trials are significantly higher than measures in baseline states. And we did not observe increased network measures in error trials.

The prefrontal cortex is at the top of the cortical perception-action cycle. WM is supported by this high-order brain area. In the Y-maze task, rats must keep reward location in free choice run “in mind” for short periods, which is used to execute subsequent alternative choices. With current multi-unit recording techniques, we can simultaneously record spikes from many neurons and track their activity changes directly. Previous studies have demonstrated that neurons in mPFC maintain WM-related information through their synaptic connections (Hains et al., 2009; Paspalas et al., 2013). The connections enable mPFC networks to efficiently process WM-related information. By studying the functional properties of mPFC networks, we may gain insight into the true nature of mPFC that underlies WM. For recent years, neurophysiologists have studied functional interactions among neurons during special cognitions (Lin et al., 2009; Stevenson et al., 2009; Hirabayashi et al., 2010; Masud and Borisyuk, 2011). In our study, well trained rats made choices at the maze choice point, on the basis of a prediction about which goal arm was the most likely rewarded one. The results show that Cd and E of mPFC networks change dynamically in correct trials. When rats started from the start point, the measures were gradually rising over time. Their peaks all appeared before the choice point. Then the measures begun to fall after rats entered the goal arms. The statistical results suggested that Cd or E in the certain durations (0.25 s pre—measure peak and 0.25 s post—measure peak) of correct trials was significantly higher than the measures in baseline states, indicating that there were real changes in Cd and E. Additionally, the significant increase pattern was not found during runs toward the choice point in error trials. mPFC network dynamics was different between correct and error trials, suggesting an association between the network states and correct behavioral choices. It is inferred that increases in mPFC network measures might have an influence as rats assessed the possible results of choices, consistent with a previous WM study (Rossetti and Carboni, 2005).

Methods for accurately estimating functional connectivity may be important for many neuroscience issues. Brain neural activity is sparse which is the point of the feature space transformation. Low-dimensional representations chosen from high dimensional data can preserve or highlight some interested features in the data (Cunningham and Yu, 2014). NMF is one of the dimensionality reduction methods. Previous studies had demonstrated that NMF could extract the dominant features in fMRI (Lohmann et al., 2007), EMG (Chiovetto et al., 2013), SPECT (Padilla et al., 2010), ECG (Ghoraani et al., 2011) database and odor profiles (Castro et al., 2013). In our study, the network analysis in NMF space presented here approaches an efficient description of mPFC network properties. The results in Figure 5A show that mPFC networks in original space are naturally sparse, in that many Granger causality values are exactly equal to zero. It is possible that some dimensions are redundant, resulting in the latent neuronal activity being confined to a low-dimensional slice of original space. NMF was used to remove dimensions that did not contribute to the raw firing rate pattern of mPFC neurons. The necessary information was correctly retained by NMF components. In Figure 5B, the average connections of mPFC networks are significantly greater than those in Figure 5A, indicating that featured mPFC activity may densely occupy the NMF space. Therefore, the functional connectivity among NMF components may offer a more accurate description of mPFC network features. The dynamics of Cd (or E) in two spaces tracking one another means that NMF components preserve the necessary information in original data. Higher Cd (or E) means that the network properties in NMF space are more apparent. The results suggested that NMF analysis faithfully captured the patterns of mPFC networks and made them more obvious. Our approach enables us to study network dynamics from a novel perspective that are not observable in original space, opening a possibility of better identifying the network characteristics.

Numerous researches have shown that mPFC is linked with WM. In the present study, food-restricted rats were rewarded if they chose the goal arm not entered in the previous free run. In this process, mPFC recalled and integrated task-relevant information to guide behavior appropriately. Previous electrophysiological evidence has suggested that lesion or inactivation of mPFC often caused increased error responding and longer retention delays in rats performing delayed spatial alternation tasks (Wang and Cai, 2006; Yoon et al., 2008; Horst and Laubach, 2009). From a molecular and cellular perspective, Kim and colleagues also have revealed that mPFC is responsible for WM by using mice with mPFC-selective PLC-β1-knockdown (Kim et al., 2014). When carried out a series of behavioral tests, they found that the mice with impaired mPFC were incapable of acquiring the task. However, rats did not depend only on the WM to perform trials during the Y-maze task. The other process might be involved. For instance, rats should inhibit immature response to the previously rewarded arm to get reward in the alternative choice run. Therefore, whether the observed mPFC activity also presumably reflects and drives many aspects of behavior (such as performance monitoring and cognitive control) deserves discussion. Previous study has suggested that the activity of a significant percentage of the mPFC neurons did represent and maintain prior task-relevant information to guide appropriate WM task behavior (Jones and Wilson, 2005). Moreover, a well-controlled study has demonstrated that elevated delay-period activity in mouse mPFC appeared to play a more important role in memory retention than inhibitory control, motor selection, and decision-making (Liu et al., 2014). These studies, together with the results reported here, support that the recorded mPFC activity closely linked with WM during Y-maze tasks.

## Conclusions

In summary, we have studied changes of mPFC spike networks during WM by means of NMF analysis. The findings point to a possible mechanism wherein the significant increases in mPFC network measures might conduct rats' subsequent choice behaviors. The anticipatory changes can be predictors of rat's future WM performance. Further results indicated that by considering the sparse feature of neural activity, the NMF method can provide a better characterization of network dynamics with less computational cost.

## Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

## Acknowledgments

Xin Tian designed the experiment; Jing Wei, Wenwen Bai and Tiaotiao Liu carried out the experiments; Jing Wei analyzed the data; Jing Wei and Xin Tian wrote the manuscript. All authors read and approved the final manuscript. This work was supported by the National Natural Science Foundation of China (91132722 and 61375113).

## Supplementary Material

The Supplementary Material for this article can be found online at: http://www.frontiersin.org/journal/10.3389/fnbeh.2015.00002/abstract

## References

Achard, S., and Bullmore, E. (2007). Efficiency and cost of economical brain functional networks. *PLoS Comput. Biol*. 3:e17. doi: 10.1371/journal.pcbi.0030017

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Anderson, K. L., Rajagovindan, R., Ghacibeh, G. A., Meador, K. J., and Ding, M. (2010). Theta oscillations mediate interaction between prefrontal cortex and medial temporal lobe in human memory. *Cereb. Cortex* 20, 1604–1612. doi: 10.1093/cercor/bhp223

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Baddeley, A. (2003). Working memory: looking back and looking forward. *Nat. Rev. Neurosci*. 4, 829–839. doi: 10.1038/nrn1201

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Baeg, E. H., Kim, Y. B., Huh, K., Mook-Jung, I., Kim, H. T., and Jung, M. W. (2003). Dynamics of population code for working memory in the prefrontal cortex. *Neuron* 40, 177–188. doi: 10.1016/S0896-6273(03)00597-X

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Benchenane, K., Peyrache, A., Khamassi, M., Tierney, P. L., Gioanni, Y., Battaglia, F. P., et al. (2010). Coherent theta oscillations and reorganization of spike timing in the hippocampal- prefrontal network upon learning. *Neuron* 66, 921–936. doi: 10.1016/j.neuron.2010.05.013

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Benjamini, Y., and Hochberg, Y. (1995). Controlling the false discovery rate: a practical and powerful approach to multiple testing. *J. R. Stat. Soc. B (Methodol.)*. 57, 289–300.

Bernhardt, B. C., Chen, Z., He, Y., Evans, A. C., and Bernasconi, N. (2011). Graph-theoretical analysis reveals disrupted small-world organization of cortical thickness correlation networks in temporal lobe epilepsy. *Cereb. Cortex* 21, 2147–2157. doi: 10.1093/cercor/bhq291

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Bressler, S. L., Richter, C. G., Chen, Y., and Ding, M. (2007). Cortical functional network organization from autoregressive modeling of local field potential oscillations. *Stat. Med*. 26, 3875–3885. doi: 10.1002/sim.2935

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Bressler, S. L., and Seth, A. K. (2011). Wiener-Granger causality: a well established methodology. *Neuroimage* 58, 323–329. doi: 10.1016/j.neuroimage.2010.02.059

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Castro, J. B., Ramanathan, A., and Chennubhotla, C. S. (2013). Categorical dimensions of human odor descriptor space revealed by non-negative matrix factorization. *PLoS ONE* 8:e73289. doi: 10.1371/journal.pone.0073289

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Chiovetto, E., Berret, B., Delis, I., Panzeri, S., and Pozzo, T. (2013). Investigating reduction of dimensionality during single-joint elbow movements: a case study on muscle synergies. *Front. Comput. Neurosci*. 7:11. doi: 10.3389/fncom.2013.00011

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Constantinidis, C., and Goldman-Rakic, P. S. (2002). Correlated discharges among putative pyramidal neurons and interneurons in the primate prefrontal cortex. *J. Neurophysiol*. 88, 3487–3497. doi: 10.1152/jn.00188.2002

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Cunningham, J. P., and Yu, B. M. (2014). Dimensionality reduction for large-scale neural recordings. *Nat. Neurosci*. 17, 1500–1509. doi: 10.1038/nn.3776

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Demirci, O., Stevens, M. C., Andreasen, N. C., Michael, A., Liu, J., White, T., et al. (2009). Investigation of relationships between fMRI brain networks in the spectral domain using ICA and Granger causality reveals distinct differences between schizophrenia patients and healthy controls. *Neuroimage* 46, 419–431. doi: 10.1016/j.neuroimage.2009.02.014

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Estrada, E., and Hatano, N. (2008). Communicability in complex networks. *Phys. Rev. E. Stat. Nonlin. Soft. Matter Phys*. 77(3 Pt 2):036111. doi: 10.1103/PhysRevE.77.036111

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Frye, R. E., Liederman, J., McGraw Fisher, J., and Wu, M. H. (2012). Laterality of temporoparietal causal connectivity during the prestimulus period correlates with phonological decoding task performance in dyslexic and typical readers. *Cereb. Cortex* 22, 1923–1934. doi: 10.1093/cercor/bhr265

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Ghoraani, B., Krishnan, S., Selvaraj, R. J., and Chauhan, V. S. (2011). T wave alternans evaluation using adaptive time-frequency signal analysis and non-negative matrix factorization. *Med. Eng. Phys*. 33, 700–711. doi: 10.1016/j.medengphy.2011.01.007

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Hains, A. B., Vu, M. A., Maciejewski, P. K., van Dyck, C. H., Gottron, M., and Arnsten, A. F. (2009). Inhibition of protein kinase C signaling protects prefrontal cortex dendritic spines and cognition from the effects of chronic stress. *Proc. Natl. Acad. Sci. U.S.A*. 106, 17957–17962. doi: 10.1073/pnas.0908563106

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Hirabayashi, T., Takeuchi, D., Tamura, K., and Miyashita, Y. (2010). Triphasic dynamics of stimulus-dependent information flow between single neurons in macaque inferior temporal cortex. *J. Neurosci*. 30, 10407–10421. doi: 10.1523/JNEUROSCI.0135-10.2010

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Horst, N. K., and Laubach, M. (2009). The role of rat dorsomedial prefrontal cortex in spatial working memory. *Neuroscience* 164, 444–456. doi: 10.1016/j.neuroscience.2009.08.004

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Hwang, K., Hallquist, M. N., and Luna, B. (2013). The development of hub architecture in the human functional brain network. *Cereb. Cortex* 23, 2380–2393. doi: 10.1093/cercor/bhs227

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Jo, Y. S., Lee, J., and Mizumori, S. J. (2013). Effects of prefrontal cortical inactivation on neural activity in the ventral tegmental area. *J. Neurosci*. 33, 8159–8171. doi: 10.1523/JNEUROSCI.0118-13.2013

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Jones, M. W., and Wilson, M. A. (2005). Theta rhythms coordinate hippocampal-prefrontal interactions in a spatial memory task. *PLoS Biol*. 3:e402. doi: 10.1371/journal.pbio.0030402

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Jonides, J., Lewis, R. L., Nee, D. E., Lustig, C. A., Berman, M. G., and Moore, K. S. (2008). The mind and brain of short-term memory. *Annu. Rev. Psychol*. 59, 193–224. doi: 10.1146/annurev.psych.59.103006.093615

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Kim, S., Putrino, D., Ghosh, S., and Brown, E. N. (2011). A Granger causality measure for point process models of ensemble neural spiking activity. *PLoS Comput. Biol*. 7:e1001110. doi: 10.1371/journal.pcbi.1001110

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Kim, S. W., Seo, M., Kim, D. S., Kang, M., Kim, Y. S., Koh, H. Y., et al. (2014). Knockdown of phospholipase C-β1 in the medial prefrontal cortex of male mice impairs working memory among multiple schizophrenia endophenotypes. *J. Psychiatry Neurosci*. 39:130285. doi: 10.1503/jpn.130285

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Lee, D. D., and Seung, H. S. (1999). Learning the parts of objects by non-negative matrix factorization. *Nature* 401, 788–791. doi: 10.1038/44565

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Lewicki, M. S. (1998). A review of methods for spike sorting: the detection and classification of neural action potentials. *Network* 9, R53–R78. doi: 10.1088/0954-898X/9/4/001

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Lin, F. H., Hara, K., Solo, V., Vangel, M., Belliveau, J. W., Stufflebeam, S. M., et al. (2009). Dynamic Granger-Geweke causality modeling with application to interictal spike propagation. *Hum. Brain Mapp*. 30, 1877–1886. doi: 10.1002/hbm.20772

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Liu, D., Gu, X., Zhu, J., Zhang, X., Han, Z., Yan, W., et al. (2014). Medial prefrontal activity during delay period contributes to learning of a working memory task. *Science* 346, 458–463. doi: 10.1126/science.1256573

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Liu, X., Liu, T. T., Bai, W. W., Yi, H., Li, S. Y., and Tian, X. (2013). Encoding of rat working memory by power of multi-channel local field potentials via sparse non-negative matrix factorization. *Neurosci. Bull*. 29, 279–286. doi: 10.1007/s12264-013-1333-z

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Logothetis, N. K., and Sheinberg, D. L. (1996). Visual object recognition. *Annu. Rev. Neurosci*. 19, 577–621. doi: 10.1146/annurev.ne.19.030196.003045

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Lohmann, G., Volz, K. G., and Ullsperger, M. (2007). Using non-negative matrix factorization for single-trial analysis of fMRI data. *Neuroimage* 37, 1148–1160. doi: 10.1016/j.neuroimage.2007.05.031

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Masud, M. S., and Borisyuk, R. (2011). Statistical technique for analysing functional connectivity of multiple spike trains. *J. Neurosci. Methods* 196, 201–219. doi: 10.1016/j.jneumeth.2011.01.003

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Nedungadi, A. G., Rangarajan, G., Jain, N., and Ding, M. (2009). Analyzing multiple spike trains with nonparametric Granger causality. *J. Comput. Neurosci*. 27, 55–64. doi: 10.1007/s10827-008-0126-2

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Paatero, P., and Tapper, U. (1994). Positive matrix factorization: a non-negative factor model with optimal utilization of error estimates of data values. *Environmetrics* 5, 111–126. doi: 10.1002/env.3170050203

Padilla, P., Górriz, J. M., Ramírez, J., Lang, E. W., Chaves, R., Segovia, F., et al. (2010). Analysis of SPECT brain images for the diagnosis of Alzheimer's disease based on NMF for feature extraction. *Neurosci. Lett*. 479, 192–196. doi: 10.1016/j.neulet.2010.05.047

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Palva, J. M., Monto, S., Kulashekhar, S., and Palva, S. (2010). Neuronal synchrony reveals working memory networks and predicts individual memory capacity. *Proc. Natl. Acad. Sci. U.S.A*. 107, 7580–7585. doi: 10.1073/pnas.0913113107

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Paspalas, C. D., Wang, M., and Arnsten, A. F. (2013). Constellation of HCN channels and cAMP regulating proteins in dendritic spines of the primate prefrontal cortex: potential substrate for working memory deficits in schizophrenia. *Cereb. Cortex* 23, 1643–1654. doi: 10.1093/cercor/bhs152

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Peyrache, A., Khamassi, M., Benchenane, K., Wiener, S. I., and Battaglia, F. P. (2009). Replay of rule-learning related neural patterns in the prefrontal cortex during sleep. *Nat. Neurosci*. 12, 919–926. doi: 10.1038/nn.2337

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Reyes-Puerta, V., Sun, J. J., Kim, S., Kilb, W., and Luhmann, H. J. (2014). Laminar and columnar structure of sensory-evoked multineuronal spike sequences in adult rat barrel cortex *in vivo*. *Cereb. Cortex*. doi: 10.1093/cercor/bhu007. [Epub ahead of print].

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Rossetti, Z. L., and Carboni, S. (2005). Noradrenaline and dopamine elevations in the rat prefrontal cortex in spatial working memory. *J. Neurosci*. 25, 2322–2329. doi: 10.1523/JNEUROSCI.3038-04.2005

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Rubinov, M., and Sporns, O. (2010). Complex network measures of brain connectivity: uses and interpretations. *Neuroimage* 52, 1059–1069. doi: 10.1016/j.neuroimage.2009.10.003

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Sakaki, M., Nga, L., and Mather, M. (2013). Amygdala functional connectivity with medial prefrontal cortex at rest predicts the positivity effect in older adults' memory. *J. Cogn. Neurosci*. 25, 1206–1224. doi: 10.1162/jocn_a_00392

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Sakurai, Y. (1999). How do cell assemblies encode information in the brain? *Neurosci. Biobehav. Rev*. 23, 785–796. doi: 10.1016/S0149-7634(99)00017-2

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Sala-Llonch, R., Peña-Gómez, C., Arenaza-Urquijo, E. M., Vidal-Piñeiro, D., Bargalló, N., Junqué, C., et al. (2012). Brain connectivity during resting state and subsequent working memory task predicts behavioural performance. *Cortex* 48, 1187–1196. doi: 10.1016/j.cortex.2011.07.006

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Schlösser, R. G., Wagner, G., and Sauer, H. (2006). Assessing the working memory network: studies with functional magnetic resonance imaging and structural equation modeling. *Neuroscience* 139, 91–103. doi: 10.1016/j.neuroscience.2005.06.037

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Schweighofer, N., Doya, K., and Lay, F. (2001). Unsupervised learning of granule cell sparse codes enhances cerebellar adaptive control. *Neuroscience* 103, 35–50. doi: 10.1016/S0306-4522(00)00548-0

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Seeley, W. W., Crawford, R. K., Zhou, J., Miller, B. L., and Greicius, M. D. (2009). Neurodegenerative diseases target large-scale human brain networks. *Neuron* 62, 42–52. doi: 10.1016/j.neuron.2009.03.024

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Seth, A. K. (2010). A MATLAB toolbox for Granger causal connectivity analysis. *J. Neurosci. Methods* 186, 262–273. doi: 10.1016/j.jneumeth.2009.11.020

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Seth, A. K., Izhikevich, E., Reeke, G. N., and Edelman, G. M. (2006). Theories and measures of consciousness: an extended framework. *Proc. Natl. Acad. Sci. U.S.A*. 103, 10799–10804. doi: 10.1073/pnas.0604347103

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Shen, X., and Meyer, F. G. (2008). Low-dimensional embedding of fMRI datasets. *Neuroimage* 41, 886–902. doi: 10.1016/j.neuroimage.2008.02.051

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Sporns, O., Tononi, G., and Kötter, R. (2005). The human connectome: a structural description of the human brain. *PLoS Comput. Biol*. 1:e42. doi: 10.1371/journal.pcbi.0010042

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Stevenson, I. H., Rebesco, J. M., Hatsopoulos, N. G., Haga, Z., Miller, L. E., and Körding, K. P. (2009). Bayesian inference of functional connectivity and network structure from spikes. *IEEE Trans. Neural Syst. Rehabil. Eng*. 17, 203–213. doi: 10.1109/TNSRE.2008.2010471

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Torres-Oviedo, G., Macpherson, J. M., and Ting, L. H. (2006). Muscle synergy organization is robust across a variety of postural perturbations. *J. Neurophysiol*. 96, 1530–1546. doi: 10.1152/jn.00810.2005

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Wang, G. W., and Cai, J. X. (2006). Disconnection of the hippocampal-prefrontal cortical circuits impairs spatial working memory performance in rats. *Behav. Brain Res*. 175, 329–336. doi: 10.1016/j.bbr.2006.09.002

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

White, T., Schmidt, M., Kim, D. I., and Calhoun, V. D. (2011). Disrupted functional brain connectivity during verbal working memory in children and adolescents with schizophrenia. *Cereb. Cortex* 21, 510–518. doi: 10.1093/cercor/bhq114

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Yan, C., Gong, G., Wang, J., Wang, D., Liu, D., Zhu, C., et al. (2011). Sex- and brain size-related small-world structural cortical networks in young adults: a DTI tractography study. *Cereb. Cortex* 21, 449–458. doi: 10.1093/cercor/bhq111

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Yoon, T., Okada, J., Jung, M. W., and Kim, J. J. (2008). Prefrontal cortex and hippocampus subserve different components of working memory in rats. *Learn. Mem*. 15, 97–105. doi: 10.1101/lm.850808

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Zhou, D., Xiao, Y., Zhang, Y., Xu, Z., and Cai, D. (2013). Causal and structural connectivity of pulse-coupled nonlinear networks. *Phys. Rev. Lett*. 111:054102. doi: 10.1103/PhysRevLett.111.054102

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Zhou, Z., Chen, Y., Ding, M., Wright, P., Lu, Z., and Liu, Y. (2009). Analyzing brain networks with PCA and conditional Granger causality. *Hum. Brain Mapp*. 30, 2197–2206. doi: 10.1002/hbm.20661

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Zhu, L., Lai, Y. C., Hoppensteadt, F. C., and He, J. (2005). Characterization of neural interaction during learning and adaptation from spike-train data. *Math. Biosci. Eng*. 2, 1–23. doi: 10.3934/mbe.2005.2.1

Pubmed Abstract | Pubmed Full Text | CrossRef Full Text | Google Scholar

Keywords: working memory (WM), spikes, functional connectivity, dimensionality reduction, non-negative matrix factorization (NMF)

Citation: Wei J, Bai W, Liu T and Tian X (2015) Functional connectivity changes during a working memory task in rat via NMF analysis. *Front. Behav. Neurosci*. **9**:2. doi: 10.3389/fnbeh.2015.00002

Received: 24 September 2014; Accepted: 02 January 2015;

Published online: 30 January 2015.

Edited by:

Daeyeol Lee, Yale University School of Medicine, USACopyright © 2015 Wei, Bai, Liu and Tian. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) or licensor are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Xin Tian, School of Biomedical Engineering, Tianjin Medical University, 22# Qixiangtai Road, Tianjin, 300070, China e-mail: tianx@tmu.edu.cn