A Mathematical Model of the Olfactory Bulb for the Selective Adaptation Mechanism in the Rodent Olfactory System

To predict the odor quality of an odorant mixture, the interaction between odorants must be taken into account. Previously, an experiment in which mice discriminated between odorant mixtures identified a selective adaptation mechanism in the olfactory system. This paper proposes an olfactory model for odorant mixtures that can account for selective adaptation in terms of neural activity. The proposed model uses the spatial activity pattern of the mitral layer obtained from model simulations to predict the perceptual similarity between odors. Measured glomerular activity patterns are used as input to the model. The neural interaction between mitral cells and granular cells is then simulated, and a dissimilarity index between odors is defined using the activity patterns of the mitral layer. An odor set composed of three odorants is used to test the ability of the model. Simulations are performed based on the odor discrimination experiment on mice. As a result, we observe that part of the neural activity in the glomerular layer is enhanced in the mitral layer, whereas another part is suppressed. We find that the dissimilarity index strongly correlates with the odor discrimination rate of mice: r = 0.88 (p = 0.019). We conclude that our model has the ability to predict the perceptual similarity of odorant mixtures. In addition, the model also accounts for selective adaptation via the odor discrimination rate, and the enhancement and inhibition in the mitral layer may be related to this selective adaptation.


Introduction
Predicting the quality of an odor composed of multiple odorant components is a challenging problem. Even if the smell of each odorant component was known, the resultant quality of an odorant mixture may differ from the linear addition of the respective odorant qualities [1].
An interesting approach to predicting odorant quality is that proposed by Haddad et al., in which an odorant quality space was derived from about 1400 kinds of odorant descriptor using principal analysis [2]. Their group also applied the developed method to predict the pleasantness of an odorant [3]. Our group proposed a model that can predict glomerular Glomerular Layer. The glomerular layer, which is composed of 1805 glomerular units, takes each measured glomerular activity pattern of the odorant components as input. The number of glomerular units is consistent with the actual number of glomeruli distributed on the olfactory bulbs of mice [21]. It then generates the glomerular activity evoked by an odorant mixture taking the respiratory cycle into account. The input of the i-th glomerular unit c i (n q ) generated by odorant n p is determined by the following procedure (cf. Fig 2), where q is the index of the odorant component.
1. Obtain a whole olfactory bulb surface image of the glomerular activity pattern corresponding to odorant n q from a web database (http://gara.bio.uci.edu/).
2. Based on the procedure used in our previous study [5], the image file of the activity pattern is divided into 1805 lattices corresponding to the glomerular units. First, each pixel of the activity pattern image is converted into a z-score corresponding to its activation level based on the color bar given by the database. A lattice filter is then adapted to divide the activity pattern into 1805 lattices.
3. Calculate the average z-score of each lattice and determine the input c i (n q ).
The output I siniff,i (t) of the i-th glomerular unit at time t is given by the following equations, which are based on a respiratory function that originated from measurement data of a rabbit glomerular layer, mitral and granular layer, and a dissimilarity evaluation module. The model takes the glomerular activity patterns of odorants composing an odor as input, and considers respiration cycles to simulate the glomerular response to odorant mixture. The neural activity in mitral and granular cells is simulated based on the models proposed in a previous study [20,21]. The dissimilarity evaluation module defines a dissimilarity index E and compares the activity patterns evoked in the mitral layer by different input odorant mixtures. and was employed in the olfactory model proposed by Li et al. [19] (cf.   [19]. In addition, P ¼ ½P 1 ðo g Þ; . . . ; P k ðo g Þ; . . . ; P 1805 ðo g Þ 2 R 1805Â1 denotes the maximum output of a glomerular unit to odor o γ , and is calculated by Eq (2). This equation was derived under the assumption that a glomerular activity pattern of an odorant mixture is the linear combination of its odorant components [10]. Note that this is a simplified formulation, because Grossman et al. [15] reported that a simple linear addition is not always applicable.
( where k is the number of odorant components in an input odor, n q is the odorant component in odorant mixture o γ , and θ is a threshold variable. The threshold function f(Á) was applied to enhance the contrast, because the 2-deoxyglucose (2-DG) method can blur the spatial response. The output P i (o γ ) of the i-th glomerular unit is sent to the mitral unit. Mitral and Granular layer. The mitral and granular layer was derived from the Li-Hopfield model [19] and the Erdi model [20]. These models, which can simulate the neural dynamics caused by the interaction between mitral and granular cells, can be written as follows: Where X ¼ ½x 1 ; . . . ; x k ; . . . ; x 1805 2 R 1805Â1 and Y ¼ ½y 1 ; . . . ; y k ; . . . ; y 1805 2 R 1805Â1 are the internal states of the mitral and granular cells, respectively, I ¼ ½I 1 ; . . . ; I k ; . . . ; I 1805 2 R 1805Â1 is the input from the glomeruli, I background is the background noise, I c is the excitatory input from the olfactory cortex, and τ x = τ y = 7 [ms] are time constants [20]. H; L; G 2 R 1805Â1805 are synapse connection matrices, where H represents the connection from granular cells to mitral cells, L represents that between mitral cells, and G represents that from mitral cells to granular cells. In addition, G Y ¼ ½g y ðy 1 Þ; . . . ; g y ðy k Þ; . . . ; g y ðy 1805 Þ 2 R 1805Â1 corresponds to the membrane potentials of mitral and granular cells, respectively. Based on the Li-Hopfield model, the membrane potentials are calculated by the following equations [19]: 8 > > > < > > > : where the threshold z is set to 1.0 based on the Erdi model [20]. As several tens of mitral cells typically receive excitatory input from the same glomerulus [22], they can be considered to form a column [23]. Thus, a mitral unit represents a mitral cell column, and each mitral unit is connected to one glomerular unit. In the same manner, the granular units represent groups of granular cells. Thus, in the proposed model, the mitral units and granular units share the same spatial distribution as the glomerular units.
To determine the synapse connection matrices, we consider each type of neuron to form a layer on a spherical surface in the olfactory bulb. First, a unit is placed on a two-dimensional α-β coordination, as shown in Fig 3. Each mitral unit is then connected to the other mitral units within a distance of z m units based on the actual connection structure of mitral cells [18]. In the same manner, each granular unit is connected to the mitral units within a distance of z g units. If this distance exceeds the extent of the two-dimensional surface, the connection is folded back to the opposite end, considering the spherical surface of the olfactory bulb shown on the right of Fig 3. The gray shadow in Fig 3 is an example of the connection range from the i-th unit located at (α i , β i ). To achieve such a connection, the following equations are used to determine the synapse connection matrices: ( where d α(i,j) and d β(i,j) denote the distance between two units along the α-axis and β-axis, respectively. H 0 ¼ ½H 0ð1;1Þ ; . . . ; H 0ðj;kÞ ; . . . ; H 0ð1805;1805Þ 2 R 1805Â1805 is the synapse connection matrix from granular units to mitral units,   (10) and (11), R i,j , R w,i , R l,i are random numbers selected from the normal distribution N(1.0,0.05), and each connection parameter is divided by the connection range parameters z m or z g so that the input strength to a unit is independent of the connection range parameters. Dissimilarity Evaluation Part. The mitral cells transfer neural activity to the olfactory cortex, and higher brain functions judge the odor quality based on individual experiences [23]. In the proposed model, the dissimilarity evaluation simulates this role by calculating a dissimilarity index based on the Pearson correlation between the spatial patterns of the mitral layer evoked by different odors. The action potential of each mitral unit responding to an odor o γ (n q ) composed of odorants n q is averaged over the respiration duration (t exhale = 400 [ms]) using the following equation to give the spatial activity pattern of the mitral layer The dissimilarity index between arbitrary odors o 1 and o 2 is then calculated by the following equation, which reflects the ease of discrimination, that is, the discrimination rate obtained from the mice experiment.
Therefore, the dissimilarity index E can be compared to the observed discrimination rate to test whether the model can predict the perceptual similarity exhibited by mice, from which the selective adaptation mechanism in the olfactory system was identified.

Simulation procedure
In our simulations, the synapse connection matrix H was adjusted using a Hebbian learning rule [20,24] based on the experimental procedure applied to mice [16], and then the dissimilarity index obtained from the model was compared to the experimentally observed discrimination rate. We consider the proposed model to correctly account for selective adaptation if it can predict the perceptual similarity of mice. This section describes the simulation procedure, including a parameter adjustment and comparison method between the simulations and experiments.
Configurations. The majority of the parameters included in the proposed model were determined based on the Li-Hopfield model [19] and the Erdi model [20], as shown in Tables C and D in S1 Appendix. The newly introduced parameters in the proposed model are the threshold θ that determines the activity pattern of the glomerular layer, and z m , z g that determine the synapse connection range of mitral units and granular units, respectively. The threshold parameter was set to θ = 0.60 to enhance the strongly activated part corresponding to a zscore greater than 2 measured by the 2-DG method. The synapse connection range parameters were set to z m = 4 and z g = 15 under the assumption that granular cells have a wider influential range than the mitral cells [25]. In addition, each component of the synapse connection matrices H, L, W was determined by Eqs (10)- (12), and 20 sets of initial connection matrices were generated to perform the following simulations. The random numbers R i,j , R w,i , R l,i used to generate different connection matrices allow us to test the robustness of prediction ability against different configurations of initial values.
Prediction of perceptual similarity.

STEP 1: Simulation of conditioning training
Based on the mice experiment [16], the rewarded odor ([IA, EB, Ci], cf. Table F in S1 Appendix) was input to the proposed model. Assuming that selective adaptation is caused by the interaction between mitral cells and granular cells, the proposed model used the following Hebbian rule [20,24] to learn the synapse connection matrix H: where the learning rates were set to η 1 = 10 −5 and η 2 = 10 −3 [19], and the adjustment was terminated when the connection parameters converged to constant values. Here, the activity pattern of the mitral layer when rewarded odor [IA, EB, Ci] was input after parameter learning is denoted as S 1 ¼ ½S . . . ; S 1805;2 2 R 1805Â1 . Activity patterns S 1 evoked by the rewarded odor and S 2 evoked by the discrimination target were substituted into Eq (14), and the dissimilarity index E was calculated. As the dissimilarity index corresponds to the discrimination rate obtained from the experiments on mice, we computed the Pearson correlation between the dissimilarity indices and the discrimination rate.

Parameter analysis.
To test the dependency of the prediction accuracies on the connection structure in the mitral and granular layers, various connection ranges z m and z g were iteratively searched. Based on a previous study [25,26], we assumed that granular cells have a wider range of influence than the mitral cells, and searched for parameters in the range 3 ≦ z m ≦ 7, 13 ≦ z g ≦ 17. We then determined the combination of z m and z g that yielded the highest correlation between the dissimilarity index E and discrimination rate of rats. We also investigated the changes in the activity patterns of the mitral layer and synapse connection parameters in H along with the Hebbian learning.  Fig 4(a) shows the transformation from glomerular activity patterns of IA, EB, and Ci into the input of the glomerular layer based on the procedure described in Section 2. Fig 4(b) compares the dissimilarity index E obtained from the simulation along with the discrimination rate of mice. Fig 4(c) shows a scatter plot of the dissimilarity index E and discrimination rate. Fig 4(b) and 4(c) confirm a strong correlation (r = 0.88, p = 0.019) between the discrimination rate of mice and the dissimilarity index E.

Results
To evaluate the prediction ability of perceptual similarity, multiple testing was performed using Bonferroni's method on the discrimination rate of mice between odor pairs. The results  [27][28][29][30][31]. The middle row shows the binarized activity patterns given by adapting Eq (2) after dividing the pixels into 1805 lattices and generalizing the activity strength into the range [0, 1]. The third row shows the output from the glomerular layer, and the bottom row shows that from the mitral layer generated by the procedure described in Section 2. (b) Comparison between discrimination rates of mice obtained from experiments and dissimilarity obtained from simulations (ζ m = 4, ζ g = 15). The figure compares the dissimilarity index E obtained from the simulation and discrimination rate for each odor. The orange bars denote simulation results and the blue bars represent the experimental results. The error bars added to the experimental results address the standard deviation in 10 mice, and the error bars added to the simulation results correspond to the standard deviation are shown in Fig 4(b), where blue lines above pairs of bars denote a significant difference (p<0.01) between the corresponding odor pairs, and in the second column of Table 1, where double asterisks represent a significant difference. Multiple testing was also adopted for the dissimilarity index E between odor pairs. The results are shown in Fig 4(b) Table 1). This result is generally consistent with the discrimination rates of mice, as shown in the second column of Table 1.
In the above simulations, the Hebbian learning rule [20,24] was used to learn the synapse connection matrix H. Fig 5(a) shows the parameters in H with respect to time, and Fig 5(b) shows the activity pattern of the mitral layer with respect to the respiratory cycle. Fig 5(a) confirms the convergence of the synapse connection parameters within five respiratory cycles. Fig  5(b) demonstrates that part of the input activity pattern has been enhanced, while the remainder is inhibited as the Hebbian learning evolves.
The parameters z m and z g , determine the number of neighbor units connected to each mitral and granular unit, and correspond to the axonal length extended from the mitral cells and granular cells in the actual olfactory bulb, respectively. To test the dependency on the connection structure, the parameters z m and z g , were iteratively varied over the range 3-7 and 13-17, respectively (giving a total of 25 parameter combinations). The resultant correlations of 20 sets of synapse connection parameters. Orange and blue lines above the bars denote a significant difference of p<0.01 between odor pairs. Orange lines represent multiple comparison results from the simulation, and blue lines represent that of experiments on mice. (c) Scatter plot between discrimination rates of mice obtained from experiments and dissimilarity index E obtained from simulations. The figure shows a scatter plot between the dissimilarity index and discrimination rate. The error bars correspond to those in (b).
doi:10.1371/journal.pone.0165230.g004 Table 1. Difference between discrimination rates (Experiment) of odor pairs and that between dissimilarity index (Simulation). Olfactory Model for Selective Adaptation Mechanism between the dissimilarity index E and discrimination rate of mice [16] are shown in Fig 6, from which the best parameter combination was determined to be z m = 4 and z g = 15.

Odor pair Experiment Simulation
Discussion Fig 4(a) shows that the output of the glomerular layer for the odor [IA, EB, Ci] is a linear addition of the activity patterns of each odorant component. This representation is consistent with the measurement data provided by Belluscio and Katz [11]. The output of the mitral layer suppressed part of the activity in the glomerular layer, which may correspond to selective adaptation. Indeed, Fig 4(c) confirms that there is a strong correlation between the dissimilarity index E and the discrimination rate of mice (r = 0.88, p = 0.019). The mice find it difficult to distinguish between the full mixture odor [IA, EB, Ci] and the subset mixture odor [IA, EB]. These results can be interpreted as the mice learning to selectively adapt and attend to certain glomeruli, which happen to be very similar for the full mixture as well as the subsets that mice have difficulty identifying.   Table 1 indicate the robustness of the prediction ability of the model against random initial values of connection matrices H, W, and L. The proposed model can therefore be considered to have the ability to predict the perceptual similarity between odors. In other words, as selective adaptation was observed in the odor discrimination rate of mice and the model can predict this discrimination rate, we can conclude that the model is sufficiently accurate to account for selective adaptation. Fig 5(a) shows that the synapse connection parameters converged after Hebbian learning, and the learning step can be terminated after five respiratory cycles. Fig 5(b) shows that part of the activity pattern input from the glomerular layer to the mitral layer was enhanced, whereas the other part was suppressed. This simulation result may correspond to the measurement data reported by Giraudet et al. [19], which indicated that one of the binary mixture components generally dominates the activity pattern on the mitral layer. Therefore, the measurement data of Giraudet et al. [19] may be a neural representation of selective adaptation in the olfactory system. Further investigation into selective adaptation should therefore measure the activity of mitral cells. In addition, as the selective adaptation emerged from Hebbian learning, our approach could be used to account for the difference in perceptual characteristics caused by odor experiences.
From the results of the parameter search described in Fig 6, the connection range parameters of mitral cells and granular cells were determined to be z m = 4 and z g = 15, which yields the highest correlation with the odor discrimination rates of mice [16]. The anatomy of the olfactory bulb, however, is not consistent with this straightforward result, because the dendrites extending from the granular cells are shorter than those from the mitral cells [32]. However, as the number of granular cells is approximately 100 times that of mitral cells [26], we assumed that the granular cells have a wider influential range than the mitral cells, and searched within the range of 3 ≦ z m ≦ 7, 13 ≦ z g ≦ 17. While this assumption was derived from a previous olfactory model proposed by Linster et al. [25], it is also possible that the cortical feedback onto granular cells gives a broader view and suppresses the response of mitral cells in an odor-selective manner. It has long been known that the olfactory cortex sends feedback to the olfactory bulb [33]. The role of this feedback was analyzed by Boyd et al. [34], who found that pyramidal cells in the olfactory cortex send odor-selective excitatory feedback to glomerular cells in the olfactory bulb, which gates the response of the mitral cells.
The proposed model and the above discussion are based on the assumption that the suppression of certain activities in the olfactory bulb is the cause of selective adaptation. If this is correct, there are at least two more mechanisms that could be involved in selective adaptation.

Intraglomerular inhibition through periglomerular cells
Interglomerular inhibition can decorrelate similar sensory inputs [35]. Selective adaptation may be a side effect of such a decorrelation.

Inhibitory feedback from pyramidal cells in the olfactory cortex to the mitral cells via glomerular cells
Selective adaptation should involve a learning process and be modulated according to the given task, because successive odor discrimination experiments on mice produced substantially improved discrimination rates for odor [IA EB] against odor [IA EB Ci] [16]. Such taskoriented learning can only be directed by higher brain functions. In addition, the cortical inhibitory feedback to the olfactory bulb helps amplify the odor-evoked inhibition [34], which is an important aspect of selective adaptation. Therefore, cortical feedback is quite possibly involved in selective adaptation and its modulation. In addition to the above scenarios, Linster et al. proposed the following mechanism: 3. Synaptic modulation in the olfactory cortex Selective adaptation may also be implemented in the olfactory cortex. Combined analysis of computer models and behavior experiments on mice suggests that synaptic adaptation and synaptic potentiation in the olfactory cortex can cause odor-specific habituation [36], which is similar to selective adaptation. (Please note that this study did not discuss odorant mixtures or measure activity patterns.) Although our simulation results suggest that selective adaptation forms solely in the olfactory bulb, the above possibility should be considered in future to clarify the relationship between selective adaptation and internal representation in the olfactory bulb and olfactory cortex. To achieve this goal, further neural activity data are required from corresponding behavior experiments.

Conclusion
We have proposed an olfactory model based on previous models [19,20], as well as some biological facts revealed in previous studies. The proposed model was used to simulate an odor discrimination experiment performed by Takiguchi et al. [16], and predicted the perceptual similarity observed in mice. Our results indicate that the proposed model is able to account for selective adaptation in the olfactory system, as observed in the odor discrimination rate of mice. We plan to investigate whether the proposed model is also applicable to humans by clarifying their selective adaptation characteristics, and then improving the proposed model to predict humans' perceptual characteristics.
Supporting Information S1 Appendix. Fig A. Odor discrimination experiment using Y-maze. Fig B. Results of odor discrimination experiment.