Gamma-beta alternation in the olfactory bulb (David, Fourcaud-Trocmé et al., 2015)

 Download zip file 
Help downloading and running models
Accession:185014
This model, a simplified olfactory bulb network with mitral and granule cells, proposes a framework for two regimes of oscillation in the olfactory bulb: 1 - a weak inhibition regime (with no granule spike) where the network oscillates in the gamma (40-90Hz) band 2 - a strong inhibition regime (with granule spikes) where the network oscillates in the beta (15-30Hz) band. Slow modulations of sensory and centrifugal inputs, phase shifted by a quarter of cycle, possibly combined with short term depression of the mitral to granule AMPA synapse, allows the network to alternate between the two regimes as observed in anesthetized animals.
Reference:
1 . David F, Courtiol E, Buonviso N, Fourcaud-Trocmé N (2015) Competing Mechanisms of Gamma and Beta Oscillations in the Olfactory Bulb Based on Multimodal Inhibition of Mitral Cells Over a Respiratory Cycle. eNeuro [PubMed]
Citations  Citation Browser
Model Information (Click on a link to find other models with that property)
Model Type: Realistic Network;
Brain Region(s)/Organism: Olfactory bulb;
Cell Type(s): Olfactory bulb main mitral GLU cell; Olfactory bulb main interneuron granule MC GABA cell;
Channel(s): I_Ks;
Gap Junctions:
Receptor(s): GabaA; AMPA;
Gene(s):
Transmitter(s):
Simulation Environment: Brian; Python;
Model Concept(s): Short-term Synaptic Plasticity; Gamma oscillations; Beta oscillations; Olfaction;
Implementer(s):
Search NeuronDB for information about:  Olfactory bulb main mitral GLU cell; Olfactory bulb main interneuron granule MC GABA cell; GabaA; AMPA; I_Ks;
# -*- coding: utf-8 -*-
"""
This script provides analysis and plot functions for LFP recordings issued from multiple simulations while varying a single parameter.
Note that a large part of the oscillation analysis (time frequency map, oscillation detection) is based on OpenElectrophy 0.2

The "main" section provides the way to plot results from a simulation output file. 

Finally, it also allow to plot few spectograms of LFPs in a simulation output file.
"""

from brian import *
from oscillation_analysis import beta_gamma_detection, spectrum_analysis
import pylab as pl
import pickle
from scipy.stats import sem, nanmean, nanstd



def multi_LFP_analysis(results,osc_th=0.06,freq_cut=40.,burn_time=0.8,distinct_figures=False):
    """
    result : list of (param,LFPsig,LFPsr)
    From a list of LFPs with repetition of a set of param:
    * detect oscillations in each LFP
    * compute some oscillation properties 
    *compute and plot stats of properties across LFPs with the same param
    """

    # Analysis of data in results
    beta_stats=[]
    gamma_stats=[]
    params=[]
    all_betas=[]
    all_gammas=[]
    all_LFPlengths=[]
    for param,LFPsig,LFPsr in results: # Check each element recorded
        params.append(param) # useful to get all simulated parameters
        
        # For that element compute beta and gamma oscillations
        list_beta,list_gamma=beta_gamma_detection(LFPsig,LFPsr/Hz,freq_cut=freq_cut,osc_th=osc_th,burn_time=burn_time/second)
        all_betas.append(list_beta)
        all_gammas.append(list_gamma)
        all_LFPlengths.append(1.*LFPsig.size/LFPsr-burn_time)

    # Compute oscillation properties from all networks with the same parameter of interest, then compute mean and std of properties
    final_list_params=unique(params)
    params=array(params)
    beta_stats=[]
    gamma_stats=[]
    for par in final_list_params:
        
        valid_par=where(params==par)[0]
        
        # Prepare to record stats from each recording
        list_beta_stats=[]
        list_gamma_stats=[]
        
        # Useful function to gather parameters of interest (amp and freq) from all oscillations in a list
        def osc_stats(list_osc):
            amps,freqs,timeinoscs=[],[],[]
            for osc in list_osc:
                amps.append(osc.amplitude_max)
                freqs.append(osc.freq_max)
                dur = (osc.time_stop-osc.time_start)/all_LFPlengths[ind] # proportion of total possible time
                timeinoscs.append(dur)
            return mean(amps),mean(freqs),sum(timeinoscs),len(amps)
        
        for ind in valid_par:
            list_beta_stats.append(osc_stats(all_betas[ind]))
            list_gamma_stats.append(osc_stats(all_gammas[ind]))
            
        list_beta_stats=array(list_beta_stats)
        list_gamma_stats=array(list_gamma_stats)
        
        # Now compute stats on osc properties for all networks using the same param (use std or sem as wanted)
        beta_stats.append([[nanmean(list_beta_stats[:,ind]),nanstd(list_beta_stats[:,ind])] for ind in range(list_beta_stats.shape[1])])
        gamma_stats.append([[nanmean(list_gamma_stats[:,ind]),nanstd(list_gamma_stats[:,ind])] for ind in range(list_gamma_stats.shape[1])])
        # beta_stats.append([[nanmean(list_beta_stats[:,ind]),sem(list_beta_stats[~isnan(list_beta_stats[:,ind]),ind])] for ind in range(list_beta_stats.shape[1])])
        # gamma_stats.append([[nanmean(list_gamma_stats[:,ind]),sem(list_gamma_stats[~isnan(list_gamma_stats[:,ind]),ind])] for ind in range(list_gamma_stats.shape[1])])
        
        
    # Tranform into arrays for simpler manipulation (shape : param, osc prop, mean/std)
    beta_stats=array(beta_stats)
    gamma_stats=array(gamma_stats)

    lw=2 # linewidth for all curves

    if distinct_figures:
        figure()
    else:
        figure(10)
    ax=subplot(411) # axis for amp vs param
    errorbar(final_list_params,gamma_stats[:,0,0],yerr=gamma_stats[:,0,1],color='k',linewidth=lw)
    errorbar(final_list_params,beta_stats[:,0,0],yerr=beta_stats[:,0,1],color='r',linewidth=lw)
    ax.set_ylabel("Amplitude")
    pl.xlim(final_list_params[0],final_list_params[-1])
    pl.ylim(ymin=0.2)

    ax=subplot(412) # axis for freq vs param
    errorbar(final_list_params,gamma_stats[:,1,0],yerr=gamma_stats[:,1,1],color='k',linewidth=lw)
    errorbar(final_list_params,beta_stats[:,1,0],yerr=beta_stats[:,1,1],color='r',linewidth=lw)
    pl.ylim(15,90)
    pl.xlim(final_list_params[0],final_list_params[-1])
    ax.set_ylabel("Frequency")

    ax=subplot(413) # total time spent in gamma and in beta
    errorbar(final_list_params,gamma_stats[:,2,0],yerr=gamma_stats[:,2,1],color='k',linewidth=lw)
    errorbar(final_list_params,beta_stats[:,2,0],yerr=beta_stats[:,2,1],color='r',linewidth=lw)
    ax.set_ylabel("Length")
    pl.xlim(final_list_params[0],final_list_params[-1])
    pl.ylim(ymin=0)

    ax=subplot(414) # number of oscillations
    errorbar(final_list_params,gamma_stats[:,3,0],yerr=gamma_stats[:,3,1],color='k',linewidth=lw)
    errorbar(final_list_params,beta_stats[:,3,0],yerr=beta_stats[:,3,1],color='r',linewidth=lw)
    ax.set_ylabel("Number")
    pl.xlim(final_list_params[0],final_list_params[-1])
    pl.ylim(ymin=-1)
    return
    
def LFP_frequency_content(results,burn_time=0.8*second,freq_min=15.,freq_max=90.,verbose=False):
    """
    Compute wavelet map (or FFT ?) of each LFP, to get an average spectral content for the whole LFP
    Then compute mean and std (sem ?) of spectral profiles across LFP with the same param
    Make a plot with all average profile on the same plot
    """

    params=[]
    all_spectrums=[]
    
    for param,LFPsig,LFPsr in results: # Check each element recorded
        params.append(param) # useful to gather all simulated params
        
        Pxx,freqs=spectrum_analysis(LFPsig,LFPsr/Hz,burn_time=burn_time/second,plot_fig=False,verbose=verbose,return_full=True)
        all_spectrums.append(Pxx)
        
    all_spectrums=array(all_spectrums)
    params=array(params)
    fig=figure()
    ax=fig.add_subplot(111)
    mask_freq=(freqs>=freq_min)&(freqs<=freq_max)
    for ind,param in enumerate(unique(params)):
        mask=params==param
        ax.plot(freqs[mask_freq],all_spectrums[mask,:][:,mask_freq].mean(axis=0),label=str(param-4.),color=pl.cm.autumn(1.*(ind+1)/(len(unique(params)))),linewidth=2)
    ax.legend()
    ax.set_xlabel("Frequency (Hz)")
    ax.set_ylabel("LFP average spectrum")

def plot_few_spectrograms(results,simul_to_plot=[0],osc_th=0.06,freq_cut=40.,burn_time=0.8*second):
    """
    Detect oscillations and plot spectrogram for all simul indices in simul_to_plot
    """
            
    for ind in simul_to_plot:
        param,LFPsig,LFPsr=results[ind] # Check each element recorded
        print "param: ",param
    
        # For that element compute beta and gamma oscillations
        list_beta,list_gamma=beta_gamma_detection(LFPsig,LFPsr/Hz,freq_cut=freq_cut,osc_th=osc_th,burn_time=burn_time/second,plot_fig=True)
        
    return

if __name__=="__main__":

    distinct_figures=False

    # Files to plot
    list_filenames=["put_filename_root_here"]
    list_filenames=["test_dict_connectproba"]
    
    # Detection osc params
    osc_th=[0.2, 0.2] # oscillation threshold
    freq_cut=40.
    burn_time=0.8*second
    
    # Spectrum params
    freq_min=15.
    freq_max=90.

    for filename in list_filenames:

        # Load data
        fid=open(filename+"_multi.out","rb")
        param_dict,results=pickle.load(fid)
        fid.close()

        print "Simulation params of plotted file: ",filename 
        for key,val in param_dict.items():
            print key,': ',val
        print
        
        tmp_results=[res[:3] for res in results]  # To keep only LFP and param
        multi_LFP_analysis(tmp_results,osc_th=osc_th,freq_cut=freq_cut,burn_time=burn_time,distinct_figures=distinct_figures)
        LFP_frequency_content(tmp_results,burn_time=burn_time,freq_min=freq_min,freq_max=freq_max)
        #~ plot_few_spectrograms(tmp_results,simul_to_plot=range(6),osc_th=osc_th,freq_cut=freq_cut,burn_time=burn_time)
    
show()