Dinamica di reti neuronali - Istituto Nazionale di Fisica ... 2.4.1 Lunghezze d’onda grandi ... un...

UNIVERSITA DEGLI STUDI DI FIRENZE

ANNO ACCADEMICO 2008-2009

Facolta di Scienze Matematiche, Fisiche e Naturali

Corso di Laurea in Fisica

Tesi di Laurea

Dinamica di reti neuronali

Candidato:

Simona Olmi

Relatore: Correlatore:

prof. Roberto Livi dott. Alessandro Torcini

28 luglio 2009

Indice

1 Elementi di neurofisiologia 1

1.1 I neuroni . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1

1.2 I segnali neuronali ed il potenziale di azione . . . . . . . . . . . . . . 3

1.3 Le sinapsi . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4

1.4 Introduzione alla dinamica neuronale . . . . . . . . . . . . . . . . . . 5

1.5 La membrana cellulare del neurone . . . . . . . . . . . . . . . . . . . 7

1.6 Canali ionici e correnti ioniche . . . . . . . . . . . . . . . . . . . . . . 9

1.7 Origine del potenziale di riposo . . . . . . . . . . . . . . . . . . . . . 11

1.7.1 Un modello cellulare semplificato . . . . . . . . . . . . . . . . 11

1.7.2 L’equazione di Nerst . . . . . . . . . . . . . . . . . . . . . . . 12

1.7.3 Effetto del sodio e pompe ioniche . . . . . . . . . . . . . . . . 13

1.7.4 Effetto della permeabilita ionica e dei meccanismi di trasporto 14

1.8 Proprieta elettriche passive della membrana . . . . . . . . . . . . . . 16

2 Modelli semplificati di neuroni e reti neurali 19

2.1 Dal modello di Hodgkin-Huxley al modello leaky integrate and fire . . 20

2.2 Neuroni leaky integrate-and-fire . . . . . . . . . . . . . . . . . . . . . 25

2.2.1 Correnti sinaptiche . . . . . . . . . . . . . . . . . . . . . . . . 28

2.2.2 Adattamento del tasso di sparo e refrattarieta . . . . . . . . . 30

I

2.3 Reti neurali . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32

2.3.1 Rete di neuroni LIF globalmente accoppiati: il modello . . . . 34

2.3.2 Mappa guidata dall’evento . . . . . . . . . . . . . . . . . . . . 37

2.3.3 Stabilita lineare . . . . . . . . . . . . . . . . . . . . . . . . . . 39

2.4 Impulso di durata finita . . . . . . . . . . . . . . . . . . . . . . . . . 41

2.4.1 Lunghezze d’onda grandi . . . . . . . . . . . . . . . . . . . . . 43

2.4.2 Lunghezze d’onda piccole . . . . . . . . . . . . . . . . . . . . . 44

2.4.3 Diagramma di fase e correzioni di taglia finita . . . . . . . . . 44

2.4.4 Altri stati collettivi della rete . . . . . . . . . . . . . . . . . . 47

2.5 Desincronizzazione in reti neurali diluite . . . . . . . . . . . . . . . . 50

Capitolo 1

Elementi di neurofisiologia

In questo capitolo introdurremo alcuni concetti basilari relativi alla neurofisiologia

del neurone. In particolare, dopo una breve descrizione della morfologia del neu-

rone e delle sinapsi, analizzeremo quei particolari segnali elettrici emessi durante la

dinamica neuronale chiamati potenziali di azione i quali rappresentano l’unita ele-

mentare di informazione scambiabile fra i neuroni. Inoltre descriveremo in dettaglio

le caratteristiche della membrana cellulare, sia a riposo che eccitata, con particolare

attenzione alle sue proprieta elettriche.

1.1 I neuroni

I neuroni sono cellule del sistema nervoso altamente specializzate e dedicate al

trasferimento, immagazzinamento ed elaborazione delle informazioni. Tale elabo-

razione avviene attraverso segnali elettrici dovuti a differenze di potenziale associate

a correnti elettriche di natura ionica (essenzialmente sono rilevanti solo gli ioni so-

dio, potassio, calcio e cloro) che attraversano la membrana cellulare del neurone.

I neuroni sono in numero molto elevato e sono connessi fra di loro in modo molto

intricato (vedi fig.1.1): ad esempio nella corteccia cerebrale si ha una densita supe-

riore a 104 neuroni/mm3. Tipicamente, pur nella varieta di tipi di cellule neuronali

(neuroni della corteccia cerebrale, motoneuroni del midollo spinale, ecc.), possiamo

riconoscere tre parti morfologicamente e funzionalmente distinte: il corpo cellulare

o soma; i dendriti e gli assoni. Il soma ha una struttura compatta che risulta ap-

prossimativamente sferica (di circa 70 µm di diametro) ed e sostanzialmente l’unita

deputata all’elaborazione dell’informazione. I dendriti sono estensioni del soma con

una struttura molto ramificata che si estendono per distanze che possono raggiun-

gere il millimetro ed hanno la funzione di raccogliere i segnali provenienti dagli altri

neuroni e trasmetterli al soma. Gli assoni sono lunghe protuberanze (anche oltre

un metro nei neuroni motori degli organismi animali superiori) che si proiettano dal

soma e la loro funzione consiste nella trasmissione del segnale generato dal soma

verso i dendriti di un’altro neurone. Usando la terminologia dei circuiti elettronici si

puo quindi dire che i dendriti rappresentano il dispositivo di ingresso (”input”), gli

assoni il dispositivo di uscita (”output”) e il soma l’unita di elaborazione delle in-

formazioni. La morfologia di un particolare tipo di neurone nonche la sua posizione

nella rete neuronale forniscono indizi sulla funzione espletata: ad esempio il livello

di arborizzazione fornisce un’idea del numero di connessioni che puo ricevere e verso

quante altre cellule neuronali invii i propri segnali (vedi fig.1.2).

Figura 1.1. Neuroni della corteccia dei mammiferi osservati al microscopio. Pos-siamo distinguere neuroni con corpi cellulari triangolari e circolari: la cellula b e

un classico esempio di cellula piramidale con corpo triangolare [10].

Figura 1.2. Forme e dimensioni dei neuroni [11].

1.2 I segnali neuronali ed il potenziale di azione

Si definisce potenziale di membrana la differenza di potenziale misurata ai capi di

due elettrodi, uno posto all’interno della cellula neuronale ed uno posto nel liquido

extracellulare circostante. Quando si parla di segnale neuronale ci si riferisce alla

variazione temporale e spaziale del potenziale di membrana. Quando il neurone e

a riposo, ossia non e in qualche modo eccitato dall’esterno (e chiariremo nel seguito

cosa si intende), il potenziale di membrana assume un valore caratteristico denom-

inato potenziale di riposo, tipicamente dell’ordine di -65mV, ossia l’interno della

cellula si trova ad un potenziale inferiore rispetto all’esterno.

I potenziali di azione sono impulsi di tensione tipici generati durante la dinamica

neuronale; essi hanno una forma pressoche stereotipata e non sono soggetti ad atten-

uazione o distorsione durante la propagazione lungo l’assone. In fig.2.3 e riportata

la forma tipica di un potenziale di azione. Si notino le seguenti caratteristiche:

• l’impulso di tensione ha una durata di circa 1-2 ms ed una ampiezza misurata

fra il minimo ed il massimo di circa 100-120 mV;

• nella prima fase dell’impulso si assiste ad una crescita veolce del potenziale di

membrana fino ad arrivare ad una fase denominata di depolarizzazione dove

il potenziale di membrana diventa positivo, cioe l’interno della cellula si trova

ad un potenziale superiore rispetto all’esterno;

• nella fase di discesa l’impulso prima di ritornare al valore di riposo passa

attraverso una fase denominata di iperpolarizzazione, tipicamente della durata

di circa 10 ms (e quindi molto piu lenta della depolarizzazione), in cui la cellula

presenta un potenziale di membrana inferiore rispetto a quello di riposo.

Il potenziale di azione, che una volta generato nella cellula neuronale, viaggia lungo

l’assone ed e trasmesso agli altri neuroni, costituisce l’unita elementare associata

alla trasmissione dei segnali neuronali. Tipicamente quindi quando ci si riferisce al

segnale emesso da un neurone si intende la sequenza temporale di questi potenziali

di azione, detta anche treno di impulsi (in inglese spike train).

Figura 1.3. Forma tipica di un potenziale di azione (o impulso) [9].

1.3 Le sinapsi

La sinapsi costituisce essenzialmente la giunzione tra due neuroni ossia la struttura

attraverso la quale le informazioni sono trasferite da una cellula nervosa all’altra.

In tale contesto si definisce neurone presinaptico il neurone ”trasmettente” i poten-

ziali di azione, a monte della sinapsi, e neurone postsinaptico il neurone ”ricevente”

i potenziali di azione, a valle della sinapsi. Con questa terminologia la sinapsi e

quindi la regione in cui l’assone del neurone presinaptico ”interagisce” con il den-

drite del neurone postsinaptico. Si definisce inoltre potenziale postsinaptico (che ha

come acronimo PPS ) la risposta in tensione del neurone postsinaptico conseguente

all’arrivo del potenziale di azione proveniente dal neurone presinaptico.

Si distinguono essenzialmente due tipi di sinapsi: la sinapsi chimica e la sinap-

si elettrica (altrimenti detta gap-junction). La sinapsi chimica, della quale e ri-

portato lo schema in fig.1.4, risulta la piu comune nel cervello dei vertebrati e si

basa sul meccanismo che andiamo a descrivere. Il potenziale di azione generato

dal neurone presinaptico, giunto all’estremita dell’assone depolarizza localmente la

membrana cellulare causando il rilascio all’interno della fessura sinaptica (ossia il

piccolo spazio tra le due membrane cellulari presinaptiche e postsinaptiche), da

parte di strutture poste sull’assone denominate vescicole sinaptiche, di particolari

sostanze chimiche denominate neurotrasmettitori. Il neurotrasmettitore, non appe-

na raggiunto il versante postsinaptico della sinapsi, e rivelato da speciali molecole

(chemorecettori) poste sulla membrana postsinaptica che provocano l’apertura (o

direttamente o tramite una catena di segnali biochimici) di specifici canali attraver-

so i quali una corrente ionica fluisce dal liquido extracellulare alla cellula. L’ingresso

di questi ioni porta a sua volta ad una variazione del valore del potenziale di mem-

brana postsinaptico. Dunque in una sinapsi chimica si ha prima la trasformazione

di un segnale elettrico in un segnale chimico sulla membrana presinaptica e poi la

successiva trasformazione sulla membrana postsinaptica di un segnale chimico in un

segnale elettrico. La sinapsi elettrica realizza invece un accoppiamento elettrico tra

due neuroni attraverso canali ionici altamente specializzati (detti gap-junctions) che

collegano la membrana presinaptica e postsinaptica. La sinapsi elettrica permette

percio un flusso di corrente diretto tra neuroni adiacenti.

Figura 1.4. Esempio di sinapsi chimica: la terminazione presinaptica libera unasostanza chimica, il neurotrasmettitore, in risposta ad una depolarizzazione [11].

1.4 Introduzione alla dinamica neuronale

Abbiamo gia accennato che l’arrivo di un potenziale di azione dal neurone presinap-

tico provoca una risposta in tensione (il potenziale postsinaptico) nel potenziale di

membrana del neurone ricevente. A tal proposito si distingue tra potenziale postsi-

naptico eccitatorio (che ha come acronimo PPSE ) e potenziale postsinaptico inibito-

rio (che ha come acronimo PPSI ) a seconda che l’effetto sia quello di aumentare

ovvero di diminuire il valore del potenziale di membrana. Analogo significato ha

la distinzione tra sinapsi eccitatoria e sinapsi inibitoria ovvero tra stimolo depo-

larizzante e stimolo iperpolarizzante. Il numero di contatti sinaptici dipende dal

tipo di neurone: ad esempio i neuroni della corteccia cerebrale (neuroni corticali)

possiedono migliaia di contatti sinaptici (da 3 · 103 a 104) con gli altri neuroni della

corteccia dei quali circa l’85% sono eccitatori ed il resto inibitori. In realta solo una

frazione dell’ordine del 5-10% risultano sinapsi realmente attive [1, 13]. Per avere

un’idea dell’ampiezza di un PPSE si veda la fig.1.5 in cui sono riportati i PPSE

registrati su una cellula neuronale ed il relativo istogramma delle altezze di picco:

come si vede il valor medio e dell’ordine di 0.5 mV. In fig.1.6 e riportata una rappre-

Figura 1.5. Registrazione di 479 PPSE sul soma di cellule piramidali della cortec-cia visuale del ratto, in presenza di attivita neuronale spontanea (fig.A). Istogram-

ma delle altezze di picco dei PPSE (fig.B) [14].

sentazione schematica della dinamica neuronale che si instaura in risposta all’arrivo

di impulsi da neuroni presinaptici:

fig.A: un neurone postsinaptico i riceve impulsi da due neuroni presinaptici j = 1,2;

ui(t) e urest rappresentano rispettivamente il potenziale di membrana e il valore

del potenziale di riposo del neurone i; la quantita ǫi1(t − t(f)1 ) rappresenta il

potenziale postsinaptico generato dall’arrivo all’istante t(f)1 di un impulso dal

neurone j = 1;

fig.B: un impulso che arriva dall’altro neurone presinaptico j = 2 ad un istante t(f)2 ,

entro un intervallo di tempo sufficientemente breve, causa un secondo poten-

ziale postsinaptico che si somma al precedente; in questo regime la risposta

del neurone postsinaptico risulta approssimativamente lineare nel senso che la

risposta e circa proporzionale agli input che riceve;

fig.C: quando ui(t) raggiunge un valore tipico, θ, denominato soglia di attivazione,

il comportamento del neurone diviene altamente non lineare: e generato un

potenziale di azione (il picco dell’impulso rappresentato da una freccia e fuori

dalla scala della figura) che ha una forma stereotipata e quindi senza legame

con gli stimoli che lo hanno prodotto; inoltre il neurone, per tutta la durata

del potenziale di azione, diviene per cosı dire ”insensibile”, ovvero refrattario,

agli stimoli che gli arrivano dagli altri neuroni.

La refrattarieta neuronale si distingue in genere in refrattarieta assoluta e refrattari-

eta relativa. La refrattarieta assoluta e quell’arco temporale (di circa 2ms) corrispon-

dente alla durata del potenziale di azione in cui e impossibile che venga generato

un altro potenziale di azione. La refrattarieta relativa, che segue temporalmente la

refrattarieta assoluta, coincide con la fase di iperpolarizzazione del neurone in cui e

”difficile” ma non impossibile che il neurone venga eccitato fino ad emettere un altro

potenziale di azione. Conseguentemente il periodo di refrattarieta fornisce un limite

inferiore al minimo intervallo temporale tra due potenziali di azione consecutivi.

Figura 1.6. Rappresentazione della dinamica neuronale [9].

1.5 La membrana cellulare del neurone

La membrana cellulare del neurone (vedi fig.1.7) e composta da molecole di lipidi e

proteine. Le molecole lipidiche sono disposte in un doppio strato con uno spessore

di circa 6nm. In questa matrice lipidica si trovano alcune molecole proteiche che

attraversano tutto lo spessore della membrana cellulare entrando quindi in contatto

sia con l’interno della cellula che con il liquido extracellulare. Tali particolari pro-

teine prendono il nome di canali proteici, ovvero canali di membrana, canali acquosi

o canali ionici (la ragione di tali denominazioni sara chiara a breve). La membrana

cellulare puo essere attraversata dalle sostanze dall’esterno all’interno o viceversa

attraverso vari meccanismi. Vi sono ad esempio alcuni tipi di molecole (alcoli o

glicerolo) che attraversano la membrana sciogliendosi nel doppio strato lipidico e

riemergendo dall’altro lato; in tal caso la facilita di penetrazione dipende dal grado

di solubilita nei lipidi. Gli ioni inorganici (sodio, potassio, calcio e cloro) che, come

abbiamo gia detto, costituiscono le correnti ioniche alla base dell’attivita elettrica

neuronale, si muovono attraverso la membrana o legandosi a particolare molecole

dette molecole di trasporto, che li veicolano, ovvero attraverso i canali proteici so-

pra citati. Si noti comunque che il meccanismo di trasporto adottato durante la

generazione di un potenziale di azione risulta quello dei canali proteici, in quanto i

flussi di ioni coinvolti (dell’ordine di 106 ioni/s) risultano ben oltre le possibilita di

azione delle molecole di trasporto. Nella struttura dei canali proteici (vedi fig.1.7)

si possono riconoscere i seguenti elementi:

• un poro centrale pieno d’acqua;

• una regione del poro che agisce da filtro di selettivita regolando il transito degli

ioni in base alle dimensioni ed alle caratteristiche chimico-fisiche;

• un sistema di porte (in inglese gates) che si aprono e si chiudono in modo

stocastico facendo oscillare il canale tra uno stato di apertura ed uno stato

di chiusura; di norma lo stato di chiusura predomina quando il potenziale di

membrana si trova al valore di riposo. 1.

1Esiste peraltro l’eccezione costituita da alcuni canali per i quali lo stato di apertura e predom-inante nella membrana a riposo: si tratta per lo piu, come vedremo nel seguito, di quei canali delcloro e del potassio che risultano responsabili proprio del valore del potenziale di riposo.

Figura 1.7. a) Struttura della membrana cellulare del neurone; b) struttura di uncanale ionico [11].

1.6 Canali ionici e correnti ioniche

Quando il canale ionico si trova in uno stato di apertura il canale e detto attivato,

altrimenti quando si trova in uno stato di chiusura e detto inattivato. E oppor-

tuno sottolineare che essendo l’apertura e chiusura dei canali un processo di natura

stocastica l’attivazione o disattivazione di un canale sta a significare solo un au-

mentata o diminuita probabilita di apertura del canale e non uno stato di apertura

o chiusura continua (vedi fig.1.8). Esistono varie modalita di attivazione, ossia di

apertura, di un canale. In particolare se l’apertura del canale puo essere regolata

dal valore del potenziale di membrana i canali si dicono voltaggio-attivati. A titolo

di esempio citiamo in questa categoria di canali il canale voltaggio-dipendente del

sodio che, come vedremo in seguito (??), e quello responsabile della depolarizzazione

della membrana che provoca la fase di salita del potenziale di azione. Riguardo alla

selettivita alla specificita ionica i canali si possono distinguere in cationici e anionici

a seconda che risultino rispettivamente permeabili agli ioni positivi o negativi. La

permeabilita di una membrana ad una specie ionica, indicata con p, e una proprieta

intrinseca della membrana che misura la ”facilita” con cui gli ioni attraversano la

membrana stessa; essa e definita in modo empirico dalla relazione [1, 17]:

J = −p∆[C] (1.1)

dove J e il flusso molare (misurato in mol/(cm2·s)) e ∆[C] rappresenta la differenza

di concentrazione ionica ai due lati della membrana (misurata in mol/cm3). P ha

le dimensioni di una velocita ed e solitamente misurata in cm/s. La permeabilita

dipende solo dal tipo e dal numero di canali ionici presenti sulla membrana. In

particolare i canali cationici possono essere specifici o non specifici qualora risultino

o meno specializzati per una particolare specie ionica (ad esempio si hanno canali

specifici del Na+, K+, Ca2+). I canali anionici risultano essenzialmente costituiti

dai canali del cloro (Cl−), che e di gran lunga il maggior permeante anionico nelle

soluzioni biologiche.

Figura 1.8. Esempio di corrente di canale. Come si vede tale corrente risulta diimpulsi di forma quasi rettangolare che possono essere posti in relazione con lo

stato di apertura e chiusura dei canali stessi [15].

La conduttanza e invece una misura dell’abilita della membrana di trasportare cor-

rente elettrica ed e misurata soltitamente in Siemens, con simbolo S, dove 1S =

1Ω−1. Poiche la corrente e trasportata dagli ioni, la conduttanza di una mem-

brana non dipendera solo dalle proprieta della membrana (cioe dalla permeabilita)

ma anche dalla concentrazione ionica all’interno ed all’esterno della cellula (ovvero il

numero dei portatori liberi di carica). Non e tuttavia possibile ricavare una relazione

matematica generale tra permeabilita e conduttanza perche tale relazione dipende

strettamente dalle modalita con cui gli ioni attraversano il canale (semplice diffu-

sione attraverso i pori pieni d’acqua ovvero modelli piu complicati che analizzano

l’interazione canale-ione permeante). In generale quindi la corrente che attraversa

un canale ionico dipendera [11, 16, 17]):

• dalla conduttanza del canale;

• dal gradiente di concentrazione tra l’interno e l’esterno della cellula che tende

a produrre un flusso dalla zona a maggior concentrazione a quella a minor

concentrazione secondo la legge empirica enunciata da Fick [17]:

Jdiff = −Dd[C]

dx(1.2)

dove Jdiff rappresenta il flusso dovuto alla diffusione (misurato in numero di

ioni/(cm2·s)), D e il coefficiente di diffusione (misurato in cm2/s) e [C] e la

concentrazione ionica (qui espressa in numero di ioni/cm3);

• dalla differenza di potenziale applicata alla membrana.

1.7 Origine del potenziale di riposo

1.7.1 Un modello cellulare semplificato

In fig.1.9 e presentato un modello cellulare semplificato [11] ma che coglie l’essenza

dell’origine del potenziale di riposo della cellula neuronale. La cellula contiene ioni

potassio (K+), sodio (Na+), cloro (Cl−) ed altri anioni (A−) di grandi dimensioni ed

e a sua volta immersa in una soluzione di sodio, potassio e cloro. Le concentrazioni

ioniche mostrate in figura sono espresse in millimoli per litro (che ha per simbolo

mM dove 1mM = 10−3mol/l). Nelle cellule neuronali sono presenti altri ioni, come

il calcio e il magnesio, ma il loro contributo al potenziale di riposo e trascurabile.

In questo modello la membrana cellulare e completamente impermeabile al sodio

ed all’anione interno (A−) e risulta invece permeabile al potassio ed al cloro. Il

potassio risulta piu concentrato dentro la cellula che fuori e tende quindi a muoversi

verso l’esterno secondo il proprio gradiente di concentrazione. D’altro canto, come

abbiamo precedentemente illustrato, la superficie interna della membrana e negativa

rispetto a quella esterna, e quindi il gradiente di potenziale tende ad attrarre il

potassio dentro la cellula. In una cellula a riposo il gradiente di concentrazione ed

il gradiente elettrico sono in equilibrio. Si definisce potenziale di equilibrio di una

specie ionica il valore del potenziale di membrana al quale non si ha alcun flusso

netto di quella specie ionica. Per il cloro il gradiente di concentrazione ed il gradiente

elettrico hanno direzione opposta rispetto al potassio.

Figura 1.9. Un modello semplificato per le distribuzioni ioniche all’interno edall’esterno della cellula. Le frecce rappresentano le direzioni dei gradienti di con-centrazione ed elettrici per il potassio (K+) e per il cloro (Cl−). Le concentrazioni

sono espresse in millimoli per litro [11].

1.7.2 L’equazione di Nerst

Il potenziale di equilibrio di ogni specie ionica e legato alle concentrazioni intracel-

lulari ed extracellulari attraverso la cosiddetta equazione di Nerst :

Eione =kT

qln

[n]e[n]i

(1.3)

dove: Eione e il potenziale di equilibrio dello ione; [n]e, [n]i sono rispettivamente

le concentrazioni extracellulari ed intracellulari; k ≃ 1.38 · 10−23J/K e la costante

di Boltzmann; T e la temperatura assoluta in Kelvin; q e la carica elettrica (in

Coulomb) della specie ionica. E possibile capire l’origine dell’equazione di Nerst

attraverso il seguente ragionamento (vedi anche fig.1.10). Dalla meccanica statistica

di Boltzmann per sistemi in equilibrio termico si ha che la probabilita p(U′

< U <

U′

+ dU) che una molecola si trovi in uno stato di energia U′

< U < U′

+ dU

risulta: p(U′

< U < U′

+ dU) ∝ exp(−U/kT ) [18]. Consideriamo adesso degli

ioni con carica positiva q in un campo elettrico statico. La loro energia nel punto

x risulta U(x) = qV (x) dove V (x) e il potenziale nel punto x. La probabilita

di trovare uno ione in una regione nell’intorno del punto x e percio proporzionale

alla quantita exp(−qV (x)/kT ). Reinterpretando la densita di probabilita come una

quantita proporzionale ad una densita ionica [n(x)], ovvero ad una concentrazione,

si ha:

p(U(x1))

p(U(x2))=

[n(x1)]

[n(x2)]= exp(−q(V (x1) − V (x2))/kT ) (1.4)

e dunque all’equilibrio termodinamico si ha che un gradiente di potenziale elettro-

statico ∆V = (V (x1) − V (x2)) genera un gradiente di densita ionica . Ma poiche

questo e un assunto su uno stato di equilibrio la relazione deve valere anche nell’al-

tro senso: ossia un gradiente di concentrazione genera un gradiente di potenziale

elettrostatico. Risolvendo la (1.4) rispetto al gradiente di potenziale si ha che all’e-

quilibrio termodinamico vale ∆V = kTq

ln [n]2[n]1

che e appunto l’equazione di Nerst.

Figura 1.10. Origine del potenziale di Nerst [9].

La quantita kT/q ha le dimensioni di una differenza di potenziale elettrico e per

una specie ionica monovalente (per la quale q = e, dove e ≃ 1.6 · 10−19C e il valore

assoluto della carica dell’elettrone) e uguale a circa 25mV a temperatura ambiente

(T ≃ 300K). Per il modello cellulare mostrato in fig.1.9 nel quale il rapporto delle

concentrazioni del cloro e del potassio e uguale e pari a 1:30 si avra a T ≃ 300K:

Ek ≃ −85mV , ECl ≃ −85mV . Nel modello cellulare semplificato, essendo gli ioni

potassio e cloro gli unici in grado di muoversi attraverso la membrana ed essendo

entrambi in equilibrio a -85mV, la cellula non manifesta alcun guadagno o perdita

netta di ioni: il potenziale di riposo della cellula semplificata coincide quindi con il

potenziale di equilibrio del potassio e del cloro.

1.7.3 Effetto del sodio e pompe ioniche

Lo ione sodio (Na+) e molto piu concentrato all’esterno che all’interno della cel-

lula (vedi fig.1.9), dunque per opporsi all’ingresso del sodio dovuto al gradiente di

concentrazione la membrana deve avere una differenza di potenziale tra l’interno e

l’esterno positiva: infatti usando le concentrazioni di fig.1.9 si ha che a T ≃ 300K

il potenziale di equilibrio del sodio risulta ENa ≃ +34mV . Dunque in una cellula

neuronale, dove il potenziale di riposo della membrana e negativo, sia il gradiente

di concentrazione sia il potenziale di membrana favoriscono l’ingresso del sodio. La

membrana cellulare e scarsamente permeabile al sodio ma il sia pur limitato in-

gresso di questo ione depolarrizza leggermente la membrana rispetto al potenziale

di equilibrio del potassio con la conseguenza che il potassio, fuori dall’equilibrio,

fluisce verso l’esterno. Per mantenere lo stato di equilibrio (ovvero le concentrazioni

ioniche costanti) a fronte di queste perdite continue esistono delle cosiddette pompe

ioniche che trasportano il sodio fuori dalla cellula ed il potassio dentro la cellula

in modo da mantenere uno stato stazionario, ovvero una sistuazione di equilibrio

dinamico. Il piu importante sistema di trasporto e costituito dalla pompa Na-K che

trasporta tre ioni sodio fuori dalla cellula per ogni due ioni potassio portati invece

all’interno: in tal caso si dice che il rapporto di accoppiamento Na:K della pompa e

di 3:2. Tale situazione e rappresentata in fig.1.11 dove le frecce tratteggiate indicano

i movimenti ionici passivi (ovvero guidati solo dal bilancio tra i gradienti elettrici e

di concentrazione) mentre le frecce continue e i circoli indicano le pompe ioniche. La

lunghezza delle frecce indica l’entita dei movimenti netti ionici: per ciascuno ione il

flusso totale e nullo all’equilibrio.

Figura 1.11. Flussi ionici passivi e pompe ioniche in una cellula in statostazionario [11].

1.7.4 Effetto della permeabilita ionica e dei meccanismi di

trasporto

Abbiamo visto nel paragrafo 1.7.2 che il potenziale di riposo della membrana e

determinato principalmente dal rapporto delle concentrazioni transmembranali del

potassio. Esiste tuttavia una dipendenza non solo dalle concentrazioni ioniche ma

anche dalla permeabilita alle specie ioniche. L’equazione di campo costante, detta

anche equazione GHK dagli autori D.E. Goldman [19], A.L. Hodgkin e B. Katz [20],

esprime tale dipendenza:

V m =kT

eln

pk[K+]e + pNa[Na+]e + pCl[Cl−]i

pk[K+]i + pNa[Na+]i + pCl[Cl−]e; (1.5)

dove Vm e il potenziale di riposo della membrana; e (valore assoluto della cari-

ca dell‘elettrone) indica la carica degli ioni monovalenti; il simbolo [I]j , dove I =

K+,Na+,Cl− e j = i,e, indica la concentrazione interna ed esterna delle specie

ioniche; pk, pNa, pCl rappresentano le permeabilita ioniche di ogni ione. Tale

equazione e basata sul presupposto che, a potenziale costante, non deve cambiare

la carica sulla membrana anche se gli ioni la attraversano in piccole quantita . Di

conseguenza le correnti trasportate dalle perdite di sodio verso l’interno, di potasso

verso l’esterno e dalle altre perdite del cloro devono dare come somma zero. Altri-

menti si detetrminerebbe un accumulo, o una perdita, costante di carica e quindi una

deriva costante del potenziale di membrana. Il nome equazione di campo costante

e legato al fatto che una delle ipotesi alla base di essa e che il campo elettrico al-

l’interno della membrana sia uniforme. Come si vede l’equazione ricorda quella di

Nerst ma considera tutte le specie ioniche presenti ”pesate” con la loro permeabilita

e corrisponde pertanto a quello che intuitivamente era possibile aspettarsi: ovvero

che il potenziale di membrana tende tanto piu al valore del potenziale di equilibrio

di una specie ionica quanto piu la sua permeabilita e maggiore di quella delle altre

specie ioniche. E possibile dimostrare [11] che non vi e una grossa dipendenza del

potenziale di riposo dal cloro e pertanto il suo contributo e spesso trascurato; in tal

caso l’equazione di campo costante viene cosı riscritta:

V m =kT

eln

[K+]e + b[Na+]e[K+]i + b[Na+]i

; (1.6)

dove b = pNa/pk. Una descrizione ancora piu accurata e fornita dall’equazione di

stato stazionario [21] che considera anche gli effetti dei processi di trasporto attivo

della pompa ionica Na − K:

V m =kT

eln

r[K+]e + b[Na+]er[K+]i + b[Na+]i

; (1.7)

dove r e il rapporto di accoppiamento del sistema di trasporto (nel caso della pompa

Na − K si ha r = 3/2). Si capisce facilmente come gli effetti della pompa si

ripercuotano sul potenziale di membrana: infatti il diverso flusso di ioni Na+ e K+

(tre ioni sodio fuori dalla cellula per ogni due ioni potassio portati invece all’interno)

causato dalla pompa genera una corrente ionica netta verso l’esterno; questo flusso

di cariche positive verso l’esterno tende a iperpolarizzare la membrana ad un valore

leggermente piu negativo (ossia avvicina il potenziale di membrana un po’ di piu al

potenziale di equilibrio del potassio) di quello che ci si aspetterebbe solo in base a

meccanismi passivi. Per avere un’idea dell’entita del contrbuto delle pompe ioniche si

confronti il valore del potenziale di riposo predetto dall’equazione di campo costante

e dall’equazione di stato stazionario relativamente ad assoni di calamaro in acqua

di mare alla temperatura di T ≃ 300K: usando l’equazione di campo costante

si trova Vm ≃ −67mV , mentre usando l’equazione di stato stazionario si trova

Vm ≃ −73mV (abbiamo usato i seguenti valori realistici: [K+]i = 400mM , [K+]e =

10mM , [Na+]i = 50mM , [Na+]e = 460mM , b = pNa/pk = 0.04, r = 1.5).

1.8 Proprieta elettriche passive della membrana

Da un punto di vista elettrico la piu semplice schematizzazione di una porzione di

membrana cellulare a riposo fa ricorso a tre elementi circuitali: una resistenza (R),

una capacita (C), ed un generatore di tensione V uguale al potenziale di riposo (vedi

fig.1.12 e fig.1.13). La presenza della resistenza si spiega con l’esistenza dei canali

Figura 1.12. Schematizzazione elettrica di una porzione di membrana cellularedel neurone. Vrest indica il potenziale di riposo [1].

Figura 1.13. Equivalente circuitale dell’intera membrana cellulare. I simboli R eC in questo caso rappresentano il paralllelo rispettivamente di tutte le resistenze ele capacita delle singole porzioni di membrana; Vrest indica il potenziale di riposo;

Iinj rappresenta una corrente iniettata dentro la cellula [1].

ionici che realizzano un contatto tra l’interno e l’esterno della cellula. La resistena

di membrana e di solito riportata come una resistenza di membrana specifica, Rm,

definita come una resistenza per unita di superficie (in unita Ω · cm2). R e quindi

ottenuta dividendo Rm per l’area della membrana considerata. Rm e determinata

principalmente dalle permeabilita a riposo del potassio e del cloro. Valori tipici di

Rm variano da circa 103 Ω ·cm2 per membrane con molti canali ionici a circa 5 x 104

Ω · cm2 per membrane con pochi canali ionici. Ma oltre a permettere il passaggio

di correnti ioniche, la membrana accumula cariche sulla superfici interna ed esterna,

ed e proprio questa separazione che detemina il potenziale di membrana. La realiz-

zazione di questa separazione di cariche conferisce alla membrana le proprieta di un

condensatore. La capacita di membrana e di solito specificata in termini di una ca-

pacita di membrana specifica, Cm, definita come una capacita per unita di superficie

(in unita F/cm2). C e quindi ottenuta moltiplicando Cm per l’area della mem-

brana considerata. Cm e tipicamente dell’ordine di 1µF/cm2. Tramite la relazione

Qm = CmV , dove Qm e la quantita di carica per unita di superficie, ed assumendo

V ≃ -65mV si ha che Qm ≃ (1µ F/cm2) x (65mV) ≃ 6.5 x 10−8 C/cm2 che equiv-

ale, dividendo per il valore assoluto della carica dell’elettrone, a circa 4 x 1011 ioni

monovalenti per cm2. Ci preme sottolineare che tale schematizzazione descrive solo

il comportamento passivo della membrana e non prende in considerazione eventuali

componenti non lineari o attivi come ad esempio conduttanze voltaggio-dipendenti

(vedi ??).

Capitolo 2

Modelli semplificati di neuroni e

reti neurali

In questo capitolo uniamo tutti gli elementi costitutivi di un neurone (dendriti,

sinapsi, conduttanze dipendenti dal potenziale, assoni) al fine di ottenere un modello

funzionale della singola cellula; con questo modello a nostra disposizione possiamo

chiederci quali operazioni puo effettuare o quanto sia efficace nel codificare infor-

mazione. Iniziamo questa operazione trascurando l’albero dei dendriti e rimpiaz-

zando la descrizione dei processi di sparo fondata sulle conduttanze (ovvero la de-

scrizione fornita dalle equazioni realistiche di Hodgkin-Huxley per le correnti di mem-

brana), con una descrizione basata sull’utilizzo di modelli formali, come il modello

leaky integrate and fire, che tratteremo. Questo modello descrive in modo estrema-

mente semplificato la dinamica sottosoglia, cioe la dinamica prima dell’insorgenza

del potenziale di azione, mentre il potenziale di azione e descritto sinteticamente

come evento impulsivo, assumendo implicitamente di aver a che fare con un im-

pulso stereotipato. In questo modo possiamo ridurre incredibilmente la complessita

del problema e caratterizzare il comportamento elettrico dei neuroni attraverso una

singola equazione differenziale ordinaria piuttosto che tramite equazioni accoppiate,

non lineari e alle derivate parziali. Tali semplificazioni ci permettono di trattare reti

formate da un gran numero di neuroni interconnessi, di simularne la dinamica e di

analizzare le soluzioni collettive, non banali, e verificarne la stabilita. La compren-

sione di un qualsiasi sistema complesso passa dunque per la scelta di un livello di

descrizione che afferra le proprieta fondamentali del sistema e trascura quelle non

essenziali per lo scopo proposto.

2.1 Dal modello di Hodgkin-Huxley al modello

leaky integrate and fire

Il modello di Hodgkin-Huxley e un modello realistico a quattro dimensioni (in cui

cioe la dinamica del neurone e descritta da quattro variabili), pensato per ripro-

durre i dati sperimentali relativi ad una particolare fibra nervosa, l’assone gigante

del calamaro [5]. Questo modello e costituito da un sistema di quattro equazioni

differenziali del primo ordine che descrivono la dinamica del potenziale di membrana

e delle tre correnti ioniche fondamentali: la corrente di sodio (INa), la corrente del

potassio (IK) e la corrente di perdita, o di dispersione, (IL dall’inglese leakage cur-

rent), dovuta principalmente al cloro (Cl−), ma che riassume anche l’effetto di altre

specie ioniche presenti, ma non descritte esplicitamente. Il meccanismo alla base di

queste correnti ioniche risiede nel fatto che le conduttanze della membrana cellulare

per il sodio (gNa) e per il potassio (gK) sono dipendenti dalla differenza di poten-

ziale applicata alla membrana, dunque la probabilita che i canali si aprano aumenta

con la depolarizzazione della membrana. L’apertura dei canali avviene tuttavia con

tempi e modalita diverse: la depolarizzazione della membrana attiva inizialmente un

aumento della conduttanza del sodio (seguita da una successiva inattivazione) e solo

Figura 2.1. Fase di salita del potenziale di azione.

con un certo ritardo temporale quella del potassio. Inoltre l’effetto sul potenziale di

membrana dell’aumento della conduttanza e diverso per il sodio e per il potassio:

• per quanto riguarda il sodio si ha un processo di retroazione positiva; infatti

una piccola depolarizzazione aumenta il numero di canali aperti e questo, a

sua volta, da luogo ad un ulteriore ingresso del sodio secondo il gradiente

elettrochimico e ad una depolarizzazione ancor piu grande;

• per quanto riguarda il potassio, si ha invece che il flusso in uscita secondo

il gradiente elettrochimico, dovuto alla depolarizzazione, porta alla ripolariz-

zazione della membrana e quindi alla conseguente riduzione della conduttanza

per il potassio fino al valore di riposo (retroazione negativa).

Il profilo tipico di un potenziale di azione e quindi spiegabile, per quanto riguarda

la fase di salita, con un improvviso grande aumento della permeabilita di membrana

al sodio (vedi figura 2.1); la conseguente corrente del sodio porta rapidamente il

potenziale di membrana a spostarsi verso il potenziale di equilibrio del sodio stesso

(ENa). La fase di caduta del potenziale di azione, ovvero la ripolarizzazione della

membrana, e dovuta ad un successivo aumento della permeabilita al potassio (vedi

figura 2.2): a causa della corrente del potassio, il potenziale di membrana si sposta

Figura 2.2. Fase di caduta del potenziale di azione.

questa volta, verso il potenziale di equilibrio del potassio (Ek).

Figura 2.3. Ricostruzione teorica dell’andamento del potenziale d’azione e relativevariazioni di conduttanza del sodio e del potassio [11].

Il modello di Hodgkin-Huxley spiega inoltre l’esistenza di un livello di soglia del

potenziale di membrana al di sopra del quale nasce l’impulso; tale fenomeno puo

essere compreso immaginando di far passare una corrente attraverso la membrana

depolarizzandola solo sino alla soglia e poi di interrompere la corrente. Essendo il

potenziale di membrana lontano dal potenziale di equilibrio del potassio, vi sara un

aumento della corrente di potassio in uscita, ma verranno attivati anche alcuni canali

di sodio, aumentando la corrente di sodio in ingresso. Alla soglia questi aumenti di

corrente sono esattamente uguali e opposti (equilibrio dinamico), e la membrana e in

una condizione simile a quella di riposo, a parte il fatto che la conduttanza del sodio

ora e instabile. La depolarizzazione oltre la soglia determina cosı un aumento di

gNa sufficiente a superare immediatamente il valore della conduttanza del potassio a

riposo, ma, non appena uno ione in piu di potassio lascia la cellula, rispetto al valore

di soglia, la depolarizzazione diminuisce e la prevalenza della corrente di potassio

determina la ripolarizzazione.

Un ultimo aspetto capace di essere spiegato con questo modello e il tempo di

refrattarieta; dopo l’emissione di un potenziale di azione abbiamo a che fare con due

periodi diversi di retrattarieta: il periodo di refrattarieta assoluto e quello relativo.

Il periodo assoluto e l’intervallo durante il quale non puo essere assolutamente emes-

so un potenziale di azione, anche se il neurone e soggetto ad uno stimolo qualsivoglia

grande; il periodo relativo e l’intervallo temporale immediatamente successivo du-

rante il quale la generazione di un nuovo potenziale di azione e inibita, ma non

impossibile. Il periodo refrattario assoluto coincide essenzialmente con l’intera du-

rata del potenziale di azione ed e causato dall’inattivazione dei canali Na+ che in

precedenza si erano aperti per depolarizzare la membrana; questi canali rimangono

inattivi finche la membrana non si ripolarizza, dopodiche si chiudono, si riattivano e

riacquistano la loro capacita di aprirsi in risposta ad uno stimolo. Il periodo refrat-

tario relativo segue temporalmente quello assoluto; i canali potassio si aprono per

terminare il potenziale di azione ripolarizzando la membrana, dunque la conduttan-

za di membrana del potassio cresce drasticamente. Gli ioni K+ che fluiscono in gran

numero fuori dalla cellula fanno sı che il potenziale di membrana si avvicini molto al

potenziale di equilibrio del potassio; questo determina una breve iperpolarizzazione

della membrana, cosiche il potenziale di membrana diventa piu negativo dell’usuale

potenziale di riposo. Finche la conduttanza del potassio non raggiunge il suo valore

di riposo, e necessario uno stimolo maggiore per raggiungere il valore di soglia e dar

vita ad una seconda depolarizzazione. Il raggiungimento del potenziale di equilibrio

di riposo determina la fine del periodo di refrattarieta relativo.

Un neurone dunque generera un potenziale di azione quando il suo potenziale di

membrana raggiungera un valore di soglia compreso tra −55 e −50 mV . Durante lo

sparo del potenziale di azione, il potenziale di membrana descresce rapidamente e

assume un valore che e iperpolarizzato rispetto al potenziale di soglia. Il meccanismo

per cui le conduttanze dipendenti dal potenziale K+ e Na+ producono dei potenziali

di azione e ben compreso e modellizzabile accuratamente; d’altra parte i modelli dei

neuroni possono essere semplificati se non sono esplicitamente inclusi i meccanismi

biofisici responsabili dei potenziali di azione. I modelli integrate-and-fire effettuano

questa semplificazione supponendo che il potenziale di azione sia generato ogni volta

che il potenziale di membrana del neurone modello raggiunge un valore di soglia Θ;

dopo l’emissione del potenziale di azione, al potenziale viene riassegnato un valore R

inferiore al potenziale di soglia, R < Θ. Tali modelli formali partono dal presupposto

che un potenziale di azione costituisca un evento stereotipato, per cui viene descritta

solo la dinamica sottosoglia di un’unica variabile (modello unidimensionale); inoltre

non riproducono in modo dettagliato i potenziali di azione ne includono l’effetto di

adattamento o il tempo di refrattarieta. L’insorgenza del potenziale di azione e cosı

caratterizzata solo dal tempo di sparo (o in inglese firing time).

2.2 Neuroni leaky integrate-and-fire

Il modello integrate-and-fire e stato proposto da Lapicque nel 1907, molto prima che

fosse compreso il meccanismo di generazione dei potenziali di azione; nonostante la

sua eta e la sua semplicita, questo modello rappresenta ancora una descrizione es-

tremamente utile dell’attivita neuronale. Trascurando una descrizione biofisica del

potenziale di azione, a questi modelli rimane il compito, piu semplice, di model-

lizzare la dinamica del potenziale di membrana sotto soglia. Questo compito puo

esser svolto con vari livelli di rigore; nella versione piu semplice, sono ignorate tutte

le conduttanze di membrana attive, compresi gli input sinaptici, e l’intera condut-

tanza di membrana e schematizzata con un singolo termine di perdita, passivo,

im = gL(V − EL), con im corrente di membrana, V potenziale di membrana, EL

potenziale di riposo di un singolo ione e gL parametro da aggiustare per ottenere la

conduttanza di membrana a riposo per uno ione. Questa versione e denominata pas-

siva o modello leaky integrate-and-fire (al quale ci riferiremo anche come LIF). Le

conduttanze neuronali sono approssimativamente costanti per piccole fluttuazioni

attorno al potenziale di membrana di riposo; il modello LIF assume che questa

costanza si mantenga nell’intero range di variazione del potenziale sotto soglia. Per

alcuni neuroni, questa si rivela un’approssimazione ragionevole, mentre per altri non

lo e; con queste approssimazioni il neurone si comporta come un circuito elettrico

costituito da una resistenza ed una capacita in parallelo (figura 2.4) ed il potenziale

di membrana si determina cosı attraverso l’equazione

cmdV

dt= −im +

Ie

A(2.1)

dove cm e la capacita per unita di superficie della membrana, A e la superficie della

membrana, im e la corrente di membrana per unita di superficie ed Ie e la corrente

di elettrodo, che usualmente non e espressa come corrente per unita di superficie;

con im = gL(V − EL) si ha

cmdV

dt= −gL(V − EL) +

Ie

A. (2.2)

E conveniente moltiplicare l’equazione (2.2) per la resistenza specifica di membrana

Figura 2.4. Circuito equivalente per il modello di neurone. Il neurone e rapp-resentato a sinistra, da un singolo insieme di superficie A con una sinapsi ed unelettrodo che inietta corrente. Nel circuito equivalente a destra, il cerchio s© indicauna conduttanza sinaptica che dipende dall’attivita di un neurone presinaptico.Viene indicata inoltre una singola conduttanza gs, anche se, in genere, ce ne pos-sono essere di tipi diversi. Il cerchio v© indica una conduttanza che dipende dal

potenziale, mentre Ie e la corrente che passa attraverso l’elettrodo. [12].

rm, data in questo caso da rm = 1/gL; questo cancella il fattore gL nel membro a

destra dell’equazione e lascia un fattore cmrm = τm nel membro di sinistra, dove τm

e la costante di tempo di membrana del neurone. L’ultimo termine del membro di

destra varia in accordo al fatto che rm/A = Rm, dove Rm e la resistenza di membrana

totale. Si arriva cosı all’equazione fondamentale per i modelli leaky integrate-and-

fire:

τmdV

dt= EL − V + RmIe. (2.3)

Per generare i potenziali di azione in questo modello occorre aggiungere una regola:

assegnare al potenziale il valore R ogniqualvolta il potenziale V raggiunga il valore di

soglia e venga sparato un potenziale di azione. L’equazione (2.3) indica che, quando

Ie = 0, il potenziale di membrana si rilassa esponenzialmente con una costante di

tempo τm a V = EL; EL e dunque il potenziale di riposo della cellula modello.

Il potenziale di membrana per il modello passivo integrate-and-fire e determinato

integrando l’equazione (2.3) e applicando la regola di soglia e di riassegnazione del

valore del potenziale per la generazione del potenziale di azione. In figura (2.5) e

mostrata la risposta di un neurone LIF ad una corrente di elettrodo che varia nel

tempo.

Figura 2.5. Modello leaky integrate and fire guidato da una corrente di elettrodoche varia temporalmente. Il grafico superiore indica il potenziale di membrana,mentre quello inferiore la corrente guida. I potenziali d’azione sono semplicementesovrapposti alla traiettoria del potenziale di membrana ogniqualvolta il potenziale

raggiunga il valore di soglia. [12].

Il tasso di sparo di un modello LIF in risposta ad una corrente iniettata costante

puo essere calcolato analiticamente. Quando Ie e indipendente dal tempo, il poten-

ziale sotto soglia V (t) puo essere facilmente calcolato risolvendo l’equazione (2.3):

V (t) = EL + RmIe + (V (0) − EL − RmIe)exp(−t/τm), (2.4)

dove V (0) e il valore di V al tempo t = 0. Questa equazione e valida per il modello

integrate-and-fire solo finche V rimane sotto soglia. Supponiamo che a t = 0 il

neurone abbia appena sparato un potenziale di azione e si abbia V (0) = R; il

potenziale di azione successivo si avra quando il potenziale di membrana raggiunge

il valore di soglia, cioe, ad un tempo t = tisi dove

V (tisi) = Θ = EL + RmIe + (R − EL − RmIe)exp(−tisi/τm). (2.5)

Risolvendo per tisi, ovvero per il tempo a cui si ha il potenziale di azione successivo,

possiamo determinare l’intervallo tra uno sparo e l’altro (in inglese interspike interval

o ISI) per una corrente Ie costante, o, equivalentemente, il suo inverso, che possiamo

chiamare il tasso di sparo del neurone,

risi =1

tisi=(

τm ln(RmIe + El − R

RmIe + EL − Θ))−1

. (2.6)

Questa espressione e valida se RmIe > Θ − EL, altrimenti risi = 0. Per valori

sufficientemente grandi di Ie, possiamo approssimare linearmente il logaritmo ln(1+

z) ≈ z per piccoli z in modo da avere

risi ≈[

RmIe + El − Θ

τm(Θ − R)

]

+

, (2.7)

che mostra che il tasso di sparo cresce linearmente con Ie per grandi Ie.

La figura (2.6A) mostra risi in funzione di Ie; il tasso di sparo, determinato come

l’inverso dell’intervallo temporale tra i primi due impulsi sparati da un neurone

corticale in vivo in risposta ad una corrente iniettata (cerchi pieni), e in accordo con i

risultati ricavati dal modello LIF, anche se il neurone reale esibisce una caratteristica

di adattamento del tasso di sparo prima di raggiungere uno stato stazionario (figura

2.6B), che il modello LIF puo riprodurre solo usando parametri diversi da quelli che

venivano usati per riprodurre gli impulsi iniziali.

2.2.1 Correnti sinaptiche

Nel caso in cui il neurone riceva degli impulsi da altri neuroni, la cosa si puo mod-

elizzare nell’ambito dei modelli LIF sostituendo alla corrente esterna Ie che compare

nell’equazione (2.3) la corrente sinaptica Isin. In particolare l’effetto prodotto da un

solo impulso ricevuto dal neurone al tempo t corrisponde ad una corrente

Isin = −Gm[V + W ]τmδ(t) (2.8)

per una sinapsi che risponda in modo infinitamente rapido [29]. Dove Gm = 1/Rm e

la conduttanza di membrana, W e il potenziale di inversione della sinapsi il cui segno

Figura 2.6. (A)Comparazione dei tassi di sparo in funzione della corrente iniet-tata, calcolati l’uno con il modello integrate and fire, l’altro attraverso le misuresu un neurone corticale in vivo. I dati puntiformi derivano da una cellula pirami-dale nella corteccia visiva primaria di un gatto. I cerchi pieni mostrano l’inversodell’ISI per i primi due impulsi sparati, mentre i cerchi vuoti mostrano il tassodi sparo nello stato quasi stazionario, dopo l’adattamento del tasso di sparo. (B)Registrazione dello sparo di un neurone corticale, soggetto ad una corrente ini-ettata costante, in cui e visibile l’adattamento. (C) Traiettoria del potenziale dimembrana e potenziali di azione per un modello integrate and fire con una corrente

aggiuntiva. [12].

determina se l’effetto dell’impulso e eccitatorio o inibitorio. Dato che il potenziale

di membrana V ha segno negativo, la sinapsi e detta eccitatoria (rispettivamente

inibitoria) se W < 0 (W > 0); il segno di V e all’origine anche del segno negativo

nell’equazione (2.8). Nel caso in cui il neurone sia connesso ad altri N neuroni

pre-sinaptici il modello LIF si puo riscrivere come

τmV = EL − V − τm(V + W )N∑

j=1

∑

m

δ(t − t(m)j ) , (2.9)

dove il tempo t(m)j rappresenta l’istante di ricezione dell’impulso m-esimo emesso dal

j-esimo neurone pre-sinaptico.

Un modello ancora piu semplificato, ma largamente usato in letteratura e che

mette in connessione la dinamica LIF con quella di oscillatori e il seguente

v = c − v −N∑

j=1

∑

m

δ(t − t(m)j ) , (2.10)

dove si e assunto che EL = 0, che il tempo unitario sia τm, che il potenziale di

membrana sia stato riscalato come v = (V − R)/(Θ − R), in modo che v ∈ [0,1].

Nell’equazione 2.10 il termine c rappresenta una corrente esterna adimensionale. In

assenza di stimoli esterni, la soluzione dell’equazione 2.10 e semplicemente v(t) =

150 175 200tempo

-60

-58

-56

-54

-52

-50

v

Θ

R

Figura 2.7. Andamento temporale del potenziale di membrana per un neuroneLIF soprasoglia in assenza di stimoli presinaptici. Le frecce indicano l’istante in

cui avviene l’emissione di un potenziale di azione.

c(1 − exp(−t)), quindi se c < 1 il suo valore asintotico resta sempre sotto soglia,

altrimenti se c > 1 il potenziale di membrana raggiunge soglia e viene riazzerato in

modo ripetitivo con periodo tisi = ln[(c − 1)/c]. In questo caso il neurone e detto

soprasoglia e la sua dinamica ricorda quella di un semplice oscillatore armonico

sovrasmorzato ed e riportata in figura 2.7

2.2.2 Adattamento del tasso di sparo e refrattarieta

L’adattamento del tasso di sparo e una caratteristica comune dei neuroni piramidali

corticali e consiste nell’allungamento degli ISI col passare del tempo in presenza

di una corrente costante iniettata nella cellula, prima di raggiungere un valore di

stato stazionario. La considerazione di questo fenomeno ci permette di mostrare

come un modello LIF possa essere modificato per incorporare una dinamica piu

complessa. Il modello LIF che abbiamo descritto nel paragrafo precedente e basato

su due approssimazioni separate: una descrizione altamente semplificata del poten-

ziale di azione e una approssimazione lineare per la corrente di membrana totale.

I dettagli del processo di generazione del potenziale di azione non sono importanti

per l’obiettivo proposto dalla nostra modellizzazione, inoltre la corrente di mem-

brana e modellizzata con dettagli sufficienti. Illustriamo il processo di adattamento

del tasso di sparo usando un modello di conduttanza che ha caratteristiche simili

a quelle delle conduttanze neuronali misurate che giocano un ruolo importante nel

produrre questo effetto. Modellizziamo l’adattamento del tasso di sparo includendo

nell’equazione fondamentale un’ulteriore corrente,

τmdV

dt= EL − V − rmgsra(V − EK) + RmIe. (2.11)

La conduttanza di adattamento del tasso di sparo gsra e stata schematizzata come

una conduttanza K+, cosiche, quando viene attivata, iperpolarizza il neurone, ral-

lentando qualsiasi impulso che puo sopraggiungere. EK e il potenziale di inversione

del potassio. Assumiamo che questa conduttanza si rilassi esponenzialmente a 0 con

una costante temporale τsra soddisfacendo l’equazione

τsradgsra

dt= −gsra. (2.12)

Ogniqualvolta un neurone genera un impulso, gsra e accresciuto di un fattore ∆gsra,

cioe gsra → gsra + ∆gsra; durante una sequenza di sparo ripetuta, la corrente si

costruisce attraverso una sequenza di passi causando l’adattamento del tasso di

sparo.

La probabilita che un neurone spari e significativamente ridotta nel breve lasso

temporale successivo alla generazione di un potenziale di azione; tale effetto refrat-

tario non e incluso nel modello base integrate-and-fire. Il modo piu semplice per

introdurre nel modello un periodo di refrattarieta assoluto e quello di aggiungere una

condizione alla regola base di superamento della soglia tale che proibisca lo sparo

per un periodo di tempo immediatamente successivo all’emissione dell’impulso. La

refrattarieta puo essere incorporata in un modo piu realistico aggiungendo una con-

duttanza simile a quella introdotta per l’adattamento del tasso di sparo; stavolta

la conduttanza dovra avere un tempo di recupero piu veloce ed un incremento piu

grande in seguito ad un potenziale di azione. Grazie ad un consistente incremento, il

potenziale del neurone si avvicina, in valore, al potenziale di equilibrio del potassio

EK in seguito ad un potenziale di azione, prevenendo in maniera temporanea un

ulteriore sparo e producendo un periodo di refrattarieta assoluto. Nel momento in

cui questa conduttanza approssima il valore 0, lo sparo sara possibile anche se in-

izialmente poco probabile, producendo un periodo di refrattarieta relativo. Quando

il recupero e completato, si ristabilisce la periodicita di sparo usuale.

2.3 Reti neurali

Una rete neurale e, formalmente, un grafo i cui nodi sono costituiti dai neuroni e

le cui connessioni sono costituite dalle sinapsi; ciascuna connessione e pesata con

la corrispondente efficienza sinaptica. Le sinapsi costituiscono cosı un grafo con

connessioni specifiche e direzionate. Ciascun nodo e caratterizzato da un’equazione

di evoluzione in cui lo stato del neurone dipende dai neuroni spazialmente connessi,

ovvero dai neuroni presinaptici; i pesi sinaptici possono essere fissati ma possono

anche evolvere nel tempo (plasticita sinaptica), in accordo con la storia dei due

nodi connessi dalla sinapsi. Se consideriamo le reti neurali come sistemi dinamici,

possiamo dare una formulazione canonica della dinamica neuronale; ciascun neurone

i e caratterizzato dal suo stato, Xi, che appartiene ad un certo insieme compatto I ∈

RM. M e il numero di variabili che caratterizzano lo stato di un neurone; assumiamo

inoltre che tutti i neuroni siano descritti dallo stesso numero di variabili. Il LIF

corrisponde al caso in cui M = 1, Xi = Vi e il potenziale di membrana del neuroni

i e I = [Vmin,Vmax]. L’evoluzione di N neuroni e data cosı da un sistema dinamico

deterministico del tipo:

dX

dt= Fγ(X,t), tempo continuo, (2.13)

o,

X(t + 1) = Fγ[X(t),t], tempo discreto. (2.14)

La variabile X = XiNi=1 rappresenta lo stato dinamico di una rete con N neuroni

al tempo t; tipicamente X ∈ M = IN , dove M e lo spazio delle fasi dell’equazione

(2.14) ed Fγ(M) ⊂ M. La mappa Fγ : M → M dipende da un set di parametri

γ ∈ RP ; il caso tipicamente considerato in seguito sara γ = (W,I(ext)), dove W e

la matrice dei pesi sinaptici, mentre I(ext) rappresenta una corrente o uno stimolo

esterno. Cosı γ e un punto nello spazio dei parametri di controllo P = N2 + N

dimensionale.

Un neurone i emette un potenziale d’azione ogni qual volta il suo stato Xi ap-

partiene ad una certa regione connessa P1 del suo spazio delle fasi, altrimenti e qui-

escente se X ∈ P0 = I\P1. Per N neuroni identici questo porta ad una “partizione

naturale” P dello spazio delle fasi M. Chiamiamo Λ = 0,1N , ω= [ωi]Ni=1 ∈ Λ,

allora P = Pωω∈Λ, dove Pω = Pω1×Pω2

× · · · × PωN. Equivalentemente, se X

∈ Pω, si ha che tutti i neuroni con ωi = 1 stanno sparando, mentre quelli con ωk = 0

sono quiescenti.

A ciascuna condizione iniziale X ∈ M possiamo associare un “raster plot” ω = ω

(t)+∞t=0 tale che X (t) ∈ Pω(t), ∀t ≥ 0; cosı ω e la sequenza di configurazioni di sparo

mostrata dalla rete neurale quando viene preparata nella condizione iniziale X (vedi

figura 2.8). Conoscere il raster plot e dunque equivalente a conoscere la lista dei

tempi di sparo di tutti i neuroni. Ad esempio, se abbiamo a che fare con una rete

omogenea di N neuroni, dove la dinamica di ciascun neurone e periodica, possiamo

dedurre, dall’osservazione del raster plot, che la dinamica della rete e asincrona;

226 226.2 226.4 226.6 226.8 227 227.2Tempo

0

0.2

0.4

0.6

0.8

1

vk

1750 1800 1850 1900 1950Tempo

0

20

40

60

80

100

Indi

ce d

i Neu

rone

Figura 2.8. (a) La figura mostra l’andamento del potenziale di membrana del k-esimo neurone in funzione del tempo, in una rete omogenea di N neuroni. (b) Nellaseconda figura e riportato l’indice del neurone che spara in funzione del tempo.[26]

dunque, se indichiamo con T l’ISI di ciascun neurone, il tasso di sparo della rete e

costante e pari a T/N .

2.3.1 Rete di neuroni LIF globalmente accoppiati: il mod-

ello

La comprensione del meccanismo di elaborazione dell’informazione nel cervello puo

essere perseguita analizzando le proprieta dinamiche dei modelli di reti neurali; pen-

sare di avvicinarsi a questo problema nella sua assoluta generalita e un obiettivo

incredibilmente ambizioso, dal momento che occorre tenere in conto (i) del ruolo

della topologia delle connessioni, (ii) della dinamica delle stesse connessioni, al fine

di rappresentare la plasticita sinaptica, (iii) della dinamica interna a ciascun mod-

ello di neurone, che puo dipendere dal numero di canali ionici ma anche da altre

variabili e parametri, (iv) dalla diversita tra i neuroni e le loro connessioni e, infine,

(v) dall’inevitabile presenza di rumore. Ad ogni modo possiamo congetturare che

almeno alcuni meccanismi di base siano robusti e possano dipendere solo da pochi

ingredienti; infatti anche i modelli semplici composti da unita identiche globalmente

accoppiate esibiscono proprieta dinamiche interessanti e, indubbiamente, non ovvie

e non totalmente comprese. Consideriamo ad esempio il cosidetto ”splay state”, che

costituisce un modo collettivo che emerge in reti di oscillatori non lineari completa-

mente accoppiati (vedi figura 2.8); questo stato e caratterizzato dal fatto che tutte

le oscillazioni abbiano la stessa forma d’onda X e che le fasi di tali oscillatori si

dispongano sul cerchio unitario in modo equispaziato. Inoltre lo stato del singolo

oscillatore xk puo essere descritto al tempo t da

xk(t) = X(t + kT/N) = A cos(ωt + 2πk/N); ω = 2π/T ; k = 1, · · · ,N (2.15)

dove N e il numero di oscillatori del sistema, T il periodo dell’oscillazione collettiva,

X la forma d’onda comune ed A l’ampiezza dell’oscillazione. La stabilita degli

stati stazionari di questi sistemi e ancora un problema dibattuto. In particolare,

nel contesto delle reti neurali, ove lo “splay state” e caratterizzato da un tasso

di sparo costante dei neuroni, si possono trovare sia dimostrazioni del fatto che

lo “splay state” sia stabile solo in presenza di accoppiamenti eccitatori [28], sia

dimostrazioni del fatto che sia possibile trovare tale stato stabile in presenza di reti

con accoppiamento completamente inibitorio [25].

Il metodo usato di solito per determinare le proprieta di stabilita di tali sem-

plici modelli e basato sull’approssimazione di campo medio; questo ci permette di

ottenere lo spettro degli autovalori associati alla matrice di stabilita nel limite ter-

modinamico N → ∞, dove N denota il numero di neuroni. La stima della stabilita

lineare dello splay state costituisce il passo preliminare verso una completa com-

prensione delle proprieta dinamiche della rete neurale; per effettuare tale stima si

puo introdurre la sezione di Poincare, vedi (??), che trasforma il sistema dinamico

originale in una mappa che connette le configurazioni dinamiche della rete neurale

corrispondenti a momenti in cui si ha l’emissione di impulsi neurali consecutivi. La

mappa viene costruita tra gli impulsi consecutivi emessi da qualsivoglia neurone e

viene combinata con un opportuno riscalamento degli indici dei neuroni in modo da

ridurre la complessita computazionale del calcolo e da ottenere espressioni analitiche

per valori grandi ma finiti di N .

Consideriamo dunque una rete di N identici neuroni LIF; la dinamica dell’i-

esimo neurone e descritta da una singola variabile, il potenziale di membrana vi(t),

che obbedisce all’equazione differenziale

vi = a − vi + gE(t) (2.16)

dove tutte le variabili e i parametri sono espressi in unita adimensionali opportuna-

mente riscalate. In accordo con la precedente equazione, il potenziale di membrana

rilassa al valore a + gE(t), ma, non appena raggiunge il valore di soglia vi = 1,

viene riscalato a vi = 0 e viene inviato, simultaneamente, un impulso a tutti gli

altri neuroni (questa procedura di riassegnazione del valore del potenziale e un mo-

do approssimato per descrivere il meccanismo di scarica presente nei neuroni reali).

Il parametro a > 1 corrisponde alla corrente di input soprasoglia, mentre g limita

la forza del campo E(t) agendo come accoppiamento efficace. Il campo E(t) e la

sovrapposizione lineare degli impulsi emessi ogniqualvolta il potenziale di membrana

di ciascun singolo neurone raggiunge il valore di soglia. Assumiamo che la forma di

un impulso emesso a t = 0 sia data da Es(t) = α2tN

e−αt , dove 1/α e l’ampiezza del-

l’impulso. Questo e equivalente ad affermare che il campo totale evolva in accordo

all’equazione

E(t) + 2αE(t) + α2E(t) =α2

N

∑

n|tn<t

δ(t − tn) , (2.17)

dove la sommatoria nel membro a destra rappresenta il termine di sorgente dovuto

agli impulsi emessi a tempi tn < t. Quindi la dinamica della rete e data da N

equazioni (2.16) piu la (2.17) che descrive l’evoluzione del campo.

2.3.2 Mappa guidata dall’evento

E conveniente trasformare le equazioni differenziali in una mappa a tempo discreto;

per fare questo occorre integrare l’equazione (2.17) tra il tempo tn e il tempo tn+1,

dove tn e il tempo immediatamente successivo all’emissione dell’n-simo impulso. La

mappa risultante e

E(n + 1) = E(n)e−ατ(n) + NQ(n)τ(n)e−ατ(n) (2.18)

Q(n + 1) = Q(n)e−ατ(n) +α2

N2, (2.19)

dove τ(n) = tn+1 − tn e l’intervallo temporale tra un impulso e il successivo e dove

abbiamo introdotto la nuova variabile Q := (αE + E)/N . L’ equazione (2.16) puo

essere integrata anche esplicitando la dipendenza temporale nota del campo E

vi(n + 1) = vi(n)e−τ(n) + a(

1 − e−τ(n))

+ gF (n) i = 1, . . . ,N . (2.20)

dove F e

F (n) =e−τ(n) − e−ατ(n)

α − 1

(

E(n) +NQ(n)

α − 1

)

− τ(n)e−ατ(n)

(α − 1)NQ(n) . (2.21)

L’intervallo temporale tra un impulso e l’altro τ e ottenuto imponendo la condizione

vm(n + 1) = 1,

τ(n) = ln

[

vm(n) − a

1 − gF (n) − a

]

. (2.22)

dove l’indice m identifica il neurone piu vicino a soglia. Le equazioni (??, 2.19)

e (2.21) possono essere scritte in una forma piu compatta esprimendo F (n) come

funzione di τ(n) per mezzo dell’equazione (2.22),

vi(n + 1) = vi(n)e−τ(n) + 1 − vm(n)e−τ(n) j = 1, . . . ,N − 1 . (2.23)

Dal momento che in una rete di neuroni indistinguibili viene preservato l’ordine dei

potenziali vi, e conveniente per prima cosa ordinarli e poi introdurre il sistema di

riferimento comovente, ovvero xj(n) = vj−n(n). In questo sistema di riferimento

l’indice del neurone piu vicino a soglia e costante e puo essere scelto sempre uguale

a 1, senza perdere di generalita; l’equazione di aggiornamento puo dunque essere

scritta come

xj−1(n + 1) = xj(n)e−τ(n) + 1 − x1(n)e−τ(n) j = 1, . . . ,N − 1 , (2.24)

con la condizione al contorno xN = 0 e dove,

τ(n) = ln

[

x1(n) − a

1 − gF (n) − a

]

. (2.25)

L’insieme di equazioni (2.18, 2.19, 2.24, 2.25, 2.21) definisce una mappa a tempo

discreto che e completamente equivalente all’insieme di equazioni differenziali di

partenza. Si ottiene cosı, che una rete di N neuroni identici puo essere descritta per

mezzo di N + 1 equazioni, due delle quali tengono conto della dinamica di E(n),

mentre le rimanenti N − 1 equazioni descrivono l’evoluzione dei neuroni (un grado

di liberta e stato perso perso andiamo da uno sparo all’altro, ovvero xN non e piu

una variabile dal momento che e stata posta, per definizione, sempre uguale a 0). Il

modello ha una dimensione finita, inoltre, in questo sistema di riferimento, lo splay

state periodico si riduce ad un punto fisso che soddisfa le seguenti condizioni:

τ(n) ≡ T

N, (2.26)

E(n) ≡ E , Q(n) ≡ Q , (2.27)

xj−1 = xje−T/N + 1 − x1e

−T/N , (2.28)

dove T e il tempo che intercorre tra due emissioni di impulso consecutive dello stesso

neurone. Attraverso un semplice calcolo si ottiene

Q =α2

N2

(

1 − e−αT/N)−1

, E = TQ(

eαT/N − 1)−1

.

La soluzione dell’equazione (2.28) coinvolge una serie geometrica, la quale, insieme

alla condizione di bordo xN = 0, porta ad un’equazione trascendentale per il periodo

T ; senza dar luogo a sorprese, il risultato e indipendente dall’ampiezza dell’impulso

α. Per semplicita riportiamo solo i termini all’ordine O(1/N) nel limite N ≫ 1,

xN−j =aT + g

T

(

1 − e−j T/N)

, (2.29)

T = ln

[

aT + g

(a − 1)T + g

]

. (2.30)

Se supponiamo di aver a che fare con neuroni disaccoppiati in regime di sparo ripeti-

tivo, cioe se prendiamo a > 1, il periodo T e ben definito nel caso eccitatorio (g > 0)

solo per g < 1 (T → 0 quando g si avvicina ad 1), mentre nel caso inibitorio (g < 0),

esiste una soluzione significativa per qualsiasi valore dell’accoppiamento (T → ∞

per g → −∞).

2.3.3 Stabilita lineare

Al fine di effettuare l’analisi di stabilita dello splay state, e necessario linearizzare le

equazioni (2.18, 2.19, 2.24) attorno al punto fisso (2.26, 2.27, 2.28),

δE(n + 1) = e−αT/NδE(n) + Te−αT/NδQ(n)

−(

αE − NQe−αT/N)

δτ(n) , (2.31)

δQ(n + 1) = e−αT/NδQ(n) − αQe−αT/Nδτ(n) , (2.32)

δxj−1(n + 1) = e−T/N [δxj(n) − δx1(n)] + e−T/N(x1 − xj)δτ(n) , (2.33)

dove abbiamo introdotto δτ(n), la cui espressione puo essere ottenuta linearizzando

le equazioni (2.25,2.21)

δτ(n) = τxδx1(n) + τEδE(n) + τQδQ(n) ; (2.34)

dove τx := ∂τ/∂x1 e dove sono state adottate definizioni analoghe per τE e τQ.

La condizione al contorno xN ≡ 0 imposta dal sistema di riferimento comovente

da δxN = 0. In pratica, il problema di stabilita e risolto calcolando lo spettro di

-1.0 -0.5 0.0 0.5 1.0Reµκ

-1.0

-0.5

0.0

0.5

1.0Imµk

Figura 2.9. Cerchio unitario e moltiplicatori esatti dello spettro di Floquet per lamappa completa per N = 20 (cerchi rossi) e N = 10 (cerchi blue) per accoppia-

mento eccitatorio. I parametri sono a = 3.0, g = 0.4, e α = 30.0.[26]

Floquet dei moltiplicatori µk, k = 1, · · · ,N + 1, associati con il problema agli au-

tovalori dell’insieme di equazioni lineari (2.31,2.32,2.33); tale calcolo deve essere, in

generale, effettuato, numericamente. Discutiamo innanzitutto il caso banale g = 0;

in questo caso si trova che µk = exp(iϕk), dove ϕk = 2πkN

, k = 1, · · · ,N − 1, e

µN = µN+1 = exp(−αT/N). Gli ultimi due esponenti interessano la dinamica del

campo di accoppiamento E(t), il cui decadimento e regolato dalla scala tempo-

rale α−1. Non appena l’accoppiamento e acceso, piccole fluttuazioni di ampiezza

∼ O(g/N) influenzano la dinamica dei neuroni e i moltiplicatori dello spettro di

Floquet assumono la forma

µk = eiϕkeT (λk+iωk)/N , ϕk =2πk

N, k = 1, . . . ,N − 1 , (2.35)

µN = eT (λN+iωN )/N , µN+1 = eT (λN+1+iωN+1)/N ,

dove λk e ωk sono la parte reale e immaginaria degli esponenti di Floquet.

In figura (2.9) mostriamo i moltiplicatori dello spettro di Floquet dello splay

state nel caso di accoppiamento eccitatorio e valori di N finiti. I moltiplicatori con

k = 1, · · · ,N − 1 sono molto vicini al cerchio unitario, mentre i due moltiplicatori

isolati µN e µN+1 giacciono molto vicini all’asse reale all’interno del cerchio unitario.

Dal momento che, gia per g/N ≈ O(10−2), i moltiplicatori del caso accoppiato

possono essere visti come una “perturbazione” piccola del caso non accoppiato (cioe

con g = 0), possiamo cercare di dar vita ad un approccio perturbativo. Prima

di effettuare questo tipo di analisi occorre sottolineare che la variabile ϕk gioca lo

stesso ruolo giocato dal numero d’onda nell’analisi di stabilita lineare di sistemi

spazialmente estesi, cosı che possiamo dire che λk caratterizza la stabilita del k-

esimo modo. Nell’analisi che segue e conveniente distinguere tra i modi caratterizzati

da ϕk ≈ 0, mod(2π) +O(1/N) e tutti gli altri modi. Si identificano dunque due

componenti spettrali che richiedono di essere trattate in modo matematicamente

diverso: la prima componente corrisponde alla condizione ‖µk − 1‖ ∼ N−1 ed e

identificata con la dicitura “lunghezze d’onda grandi” (in inglese long wavelenghts

o LWs); la seconda componente corrisponde a ‖µk − 1‖ ∼ O(1) ed e identificata con

“lunghezze d’onda brevi” (in inglese short wavelenghts o SWs).

2.4 Impulso di durata finita

Analizziamo adesso il problema di stabilita dello splay state per reti soggette ad

impulsi con α finito (cioe indipendente dalla dimensione del sistema) per grandi

valori di N . Considerando i termini all’ordine 1/N , la mappa ad evento guidato

(vedi equazioni (2.18, 2.19, 2.24) ) si semplifica nel seguente insieme di N + 1

equazioni:

E(n + 1) = (1 − ατ)E(n) + NQ(n)τ , (2.36)

Q(n + 1) = (1 − ατ)Q(n) +α2

N2, (2.37)

xj−1(n + 1) = (1 − τ)xj(n) + 1 − x1(n) + τ , (2.38)

dove j = 1, . . . ,N−1 mentre τx1 e stato approssimato con τ poiche 1−x1 ≈ O(1/N).

L’espressione per l’intervallo di intersparo (2.25) si semplifica cosı in

τ(n) =1 − x1(n)

a − 1 + gE(n). (2.39)

La soluzione periodica per il campo che tiene conto di tutti gli impulsi diventa

E = T−1 e Q = α/NT , mentre xj e il periodo T sono ancora dati dall’equazioni

(2.29, 2.30). Lo spettro degli autovalori di Floquet µk, k = 1, · · · ,N + 1 puo essere

ottenuto linearizzando le equazioni (2.36, 2.37, 2.38) nell’intorno di questa soluzione

periodica e assumendo che ciascuna perturbazione cresca come µnk

µkδE = (1 − αT/N)δE + TδQ , (2.40)

µkδQ = (1 − αT/N)δQ − α2

NTδτ , (2.41)

µkδxj−1 = (1 − T/N)δxj − δx1 + (1 − xj)δτ , (2.42)

Un’espressione esplicita per δτ si puo ottenere valutando la derivata dell’equazione

(2.39)

δτ = −T 2

NδE − T

gδx1 .

Sostituendo quanto trovato nelle equazioni (2.40,2.41,2.42), vediamo che il problema

agli autovalori consiste nel trovare le N + 1 radici µk del polinomio associato; una

semplificazione parziale del problema puo essere ottenuta estraendo dalle equazioni

(2.41) e (2.42) la dipendenza di δτ dalla perturbazione del potenziale del neurone

piu vicino a soglia δx1

δτ = Kδx1 (2.43)

dove

K =

[

−(a − 1 + g/T ) +α2gT

N2(µk − 1 + αT/N)2

]−1

. (2.44)

Infine, sostituendo questa espressione nell’equazione (2.40), otteniamo un insieme

chiuso di equazioni per le perturbazioni del potenziale di membrana,

µkδxj−1 = (1 − T/N)δxj + (K − 1)δx1 − Kxjδx1 . (2.45)

Dopo aver imposto la condizione di bordo δxN = 0, l’equazione (2.45) si riduce alla

seguente equazione agli autovalori:

µN−1k eT = K(a + g/T )

1 − µN−1k

1 − µk− [K(a − 1 + g/T ) + 1]

1 − µN−1k eT

1 − µkeT/N, (2.46)

dove K e una funzione di µk. Al fine di risolvere questa equazione analiticamente e

necessario effettuare la distinzione tra lunghezze d’onda grandi e piccole.

2.4.1 Lunghezze d’onda grandi

Consideriamo i modi per cui ‖µk−1‖ ∼ N−1 (o, equivalentemente, ϕk ≈ 0, mod(2π)

+O(1/N)). Al fine di semplificare la notazione definiamo

Λk :=N

Tln µk .

All’ordine principale, K vale

K =

[

−(a − 1 + g/T ) +α2g

T (Λk + α)2

]−1

.

Sostituendo questa espressione nell’equazione agli autovalori (2.46) e rimuovendo i

termini in 1/N , si ottiene

aT + g

α2g

(

1 − e−ΛkT)

(Λk + α)2(Λk + 1) = Λk

(

eT − e−ΛkT)

, ||Λk|| 6= 0 . (2.47)

Questa equazione coincide con quella derivata attraverso l’analisi di campo medio

[28], eccetto per il fattore ||Λk|| = 0, che corrisponde all’esponente di Floquet zero

dell’evoluzione temporale continua (2.16) e scompare nella dinamica a tempo discre-

to. Occorre sottolineare ancora che, nonostante l’equazione (2.47) abbia N+1 radici,

questa da un’opportuna approssimazione solo per quelli autovalori che soddisfano la

relazione ||µk − 1|| ∼ N−1.

2.4.2 Lunghezze d’onda piccole

La seconda componente dello spettro (2.35) e ottenuta per ||µk − 1|| ∼ O(1). In

questo caso K ha la semplice forma

K = −(a − 1 + g/T )−1 .

Introducendo nell’equazione (2.46) la forma esplicita della soluzione periodica (2.30),

si ottiene una notevole semplificazione per lo spettro, ovvero

µNk = e−T a + g/T

a − 1 + g/T= 1 , (2.48)

i.e., questo coincide con uno spettro di Floquet completamente degenere,

ωk ≡ 0 , λk ≡ 0 (2.49)

Notiamo che questa approssimazione vale solo per quelli autovalori per cui ||µk−1|| ∼

O(1), ovvero vale per la maggior parte dello spettro, eccetto quelli autovalori che

giacciono vicino al punto (1,0), dove il cerchio unitario intercetta l’asse reale (vedi

fig. 2.9).

Per quanto riguarda gli autovalori isolati µN e µN+1, si puo facilmente intuire

che questi possano essere sempre contenuti all’interno del cerchio unitario, vicino

all’asse reale, per N finiti e ampiezze dell’impulso finite, e che si avvicinino al punto

(1,0) da sinistra nel limite N → ∞; dunque questi possono al piu contribuire alla

stabilita marginale della dinamica.

2.4.3 Diagramma di fase e correzioni di taglia finita

In accordo con i precedenti risultati possiamo dire che, nel limite N → ∞, il mani-

festarsi dell’instabilita e determinato dagli esponenti di Floquet associati con le LW;

inoltre si ottiene che lo splay state e sempre instabile per accoppiamento inibitorio

0 0.2 0.4 0.6 0.8 1g

0

50

100

150

200

α

INSTABILE

STABILE

Figura 2.10. Diagramma di fase per la stabilita dello splay state in una reteneurale con accoppiamento eccitatorio che agisce attraverso un impulso di duratafinita. La linea continua che separa la regione stabile da quella instabile nel piano(g,α) e stata ricavata dalla formula analitica dello spettro di Floquet (2.47) cona = 1.3 . Occorre sottolineare che, in questo contesto, stabile si riferisce ad N finiti,

dal momento che per sistemi infiniti la stabilita diventera marginale.[26]

(g < 0). Nel caso di accoppiamento eccitatorio l’analisi di campo medio predice

stabilita per lo splay state per α ≤ αc(g,a), dove αc(g,a) e la linea critica che separa

la regione stabile da quella instabile (vedi figura 2.10). Questa linea corrisponde

ad una biforcazione di Hopf che da vita ad un nuovo comportamento periodico col-

lettivo chiamato sincronizzazione parziale . Inoltre si ha che αc diverge a +∞ per

g → 1 ( per g > 1 non puo essere sostenuto un regime stazionario dal momento che

l’evoluzione accelera in maniera uniforme ). Nel limite opposto di accoppiamento

che si annulla (g → 0), la biforcazione sopravvive; dall’equazione (2.47) si puo vedere

inoltre che la frequenza ωc della soluzione di biforcazione converge a ωc = 2π/T e

αc = −1 +√

1 + ω2c .

Includendo il ruolo delle SW possiamo concludere che, nel limite N → ∞, lo

stato stazionario puo essere al piu marginalmente stabile per α ≤ αc(g,a). La de-

generazione perfetta dell’esponente zero di Floquet associata con le lunghezze d’onda

brevi toglie i dubbi sulle proprieta di stabilita efficaci di reti grandi ma finite, dal mo-

mento che tali modi sono marginalmente stabili. Risolvendo inoltre numericamente

le equazioni agli autovalori (2.31,2.32,2.33) al variare delle dimensioni del sistema

100

101

102

103

k10

-6

10-5

10-4

10-3

10-2

|λk|(a)

-600 -400 -200 0 200 400 600ϕkN

-0.005

-0.004

-0.003

-0.002

-0.001

0λk (b)

Figura 2.11. (a) Grafico log-log dei valori assoluti degli esponenti di Floquet λk,ordinati dal piu grande al piu piccolo come funzione dell’indice k = 1,...,N perN = 100,200,400. La linea tratteggiata ha una pendenza di -2. (b) Esponentedi Floquet come funzione della fase riscalata ϕN , per N = 100 (cerchi neri) andN = 200 (cerchi rossi). In entrambe le figure i valori dei parametri sono a = 3.0,

g = 0.4, e α = 30.0.[26]

e al variare di g e α si ottiene che lo splay state e strettamente stabile in reticoli

finiti e che l’esponente di Floquet massimo tende a zero da sotto con un andamento

del tipo 1/N2 (vedi figura 2.11). Questo implica che una teoria delle perturbazioni

fatta fino all’ordine 1/N non puo tener conto delle instabilita del modello originale.

Questo e confermato dalla figura 2.12 dove gli spettri di Floquet ottenuti da ap-

prossimazioni al primo e al secondo ordine danno luogo ad uno splay state instabile,

anche se la soluzione numerica del problema di stabilita indica che il modello ad N

finito e stabile.

Riassumendo, la stabilita dello splay state puo essere indagata riducendo un

modello di equazioni differenziali globalmente accoppiate a delle opportune mappe

guidate dall’evento, che mettono in relazione le configurazioni interne della rete cor-

rispondenti a due emissioni consecutive del potenziale di azione. L’analisi analitica

dello Jacobiano nel limite di grandi N rivela che lo spettro degli autovalori e costitu-

ito da due componenti: (i) modi propri di grandi lunghezze d’onda; (ii) modi propri

di brevi lunghezze d’onda. Un approccio di campo medio e capace di riprodurre la

parte dello spettro associata ai modi di grandi lunghezze d’onda col rischio di non

catturare l’instabilita che si manifesta a lunghezze d’onda piccole.

-π -π/2 0 π/2 π ϕk

-8×10-4

-4×10-4

0

4×10-4

8×10-4

λk

Figura 2.12. Esponenti di Floquet λk(ϕ) come funzione della fase ϕk per impulsidi durata finita α = 30.0 e numero di neuroni finito N = 500 nel caso di accop-piamento eccitatorio g = 0.4. Le cerchi pieni rappresentano il risultato esatto perN finito, mentre i quadrati rossi quadrati e i triangoli blue vuoti si riferiscono arisultati approssimati corretti rispettivamente fino al primo e al secondo ordine in

1/N . Il parametro qui e a = 3.[26]

2.4.4 Altri stati collettivi della rete

Discutiamo brevemente adesso della soluzione collettiva che nasce al di sopra della

linea αc(g,a) (vedi figura 2.10), detta sincronizzazione parziale. In questo regime il

campo medio esibisce una dinamica periodica che si manifesta in assenza di qualsiasi

sincronizzazione tra i singoli neuroni, i quali si comportano in modo quasi periodico.

La sincronizzazione parziale deriva dalla destabilizzazione del regime splay state,

che e invece caratterizzato da un campo medio costante e da un comportamento

periodico dei singoli neuroni, le cui fasi sono equispaziate.

In una rete di neuroni LIF, in cui l’accoppiamento e modellizzato da una funzione

α, come in (2.17), troviamo che, per accoppiamento eccitatorio, le reti evolvono ver-

so uno stato asincrono (splay state) se le costanti di accoppiamento temporale sono

sufficientemente basse (piccoli α). Quando α aumenta, la rete inizia a sincroniz-

zarsi; per costanti di tempo finite la rete non raggiunge sincronia completa anche

se e costituita da oscillatori identici e il sistema e completamente privo di rumore.

Solo nel limite di α → ∞ e di durata dell’impulso che tende a zero, lo stato diventa

completamente sincronizzato. Dunque se attraversiamo il diagramma di fase nel

piano (g,α) mantenendo g = costante, passiamo, al crescere di α dallo splay state

alla sincronizzazione parziale e solo nel limite α → ∞ otteniamo una sicronizzazione

completa. Se la rete e composta da piu di due unita, i neuroni sparano quasi peri-

odicamente nello stato parzialmente sincronizzato; in sistemi di grandi dimensioni

il tasso di sparo medio della rete in questo stato parzialmente sincronizzato varia

periodicamente, anche se le singole unita sono quasiperiodiche.

Con accoppiamento inibitorio la rete sincronizza completamente per α piccoli;

per α grandi la rete si rompe in due o piu cluster completamente sincronizzati. Il nu-

mero medio di cluster che si forma cresce al crescere di α, sebbene il numero esatto di

cluster formati dipenda dalle condizioni iniziali. Per una rete con un numero grande

ma finito di oscillatori il numero di cluster formati tendera ad essere pari al numero

di neuroni, dunque, lo stato finale sara indistinguibile dallo splay state. Cosı con

accoppiamento inibitorio, il sistema passa da uno stato completamente sincronizza-

to per piccoli α, ad uno stato sostanzialmente asincrono per grandi α, attraverso

la formazione di un numero crescente di cluster completamente sincronizzati; nel

limite α → ∞ si ritrova lo splay state. La transizione da un sistema sincrono ad

uno asincrono e accompagnata dall’allargamento di un cluster piuttosto che dalla

rottura in piu cluster.

In figura (2.13) mostriamo alcuni risultati di simulazioni numeriche effettuate

nel regime di sincronizzazione parziale. Analizziamo in particolare la figura (b).

Se indichiamo con tn[i] il tempo a cui l’i-esimo neurone emette un impulso, e con

T (tn[i]) l’ultimo ISI dello stesso neurone, possiamo ricostriure tutti i tempi di sparo

passati iterando semplicemente la formula ricorsiva

tn−1[i] = tn[i] − T (tn[i]). (2.50)

Geometricamente, implementare questa formula corrisponde ad eseguire la seguente

procedura (fig. 2.13(b)): dato il punto (tn,T (tn)) ci muoviamo verso il basso lungo

la retta di pendenza 1 fino a raggiungere l’asse t al tempo tn−1; a questo punto ci

spostiamo verticalmente, verso l’alto, fino ad incontrare la curva T . Considerando

che T (t) e una funzione periodica definita positiva, l’equazione (2.50) non e nient’al-

tro che l’equazione di evoluzione di un oscillatore periodico forzato; ci aspettiamo

dunque l’insorgere di fenomeni di bloccaggio ( in inglese locking ). D’altra parte le

simulazioni numeriche non presentano alcun fenomeno di bloccaggio, anche quando

il rapporto tra le frequenze caratteristiche del termine forzante e del campo forzato

e razionale. Questa caratteristica rimane inspiegata.

In figura (2.13(c)) e riportato l’andamento del “flusso istantaneo” in funzione del

tempo; tale flusso e definito da

π(t) =1/N

tm+1 − tm, (2.51)

dove (tm+1−tm) ≈ (1/N) e l’intervallo tra due spari consecutivi della rete; il membro

a sinistra di questa equazione non presenta il pedice m, perche nel limite N → ∞

la variabile temporale diventa nuovamente continua.

Un’altra peculiarita dello stato di sincronizzazione parziale e il fatto che gli ISI

dei singoli neuroni differiscano dal periodo della dinamica macroscopica; in altre

parole i tempi tra un’emissione e la successiva di un impulso, da parte del singolo

neurone, sono sempre piu piccoli del periodo di oscillazione del campo forzante.

Possiamo spiegare il fenomeno, in maniera qualitativa, osservando che, da un lato,

i neuroni tendono a raggrupparsi, dall’altro, i neuroni che si trovano sul fronte del

raggruppamento tendono a scappare, mentre quelli che raggiungono il gruppo dalle

retrovie vi si uniscono. Questo giustifica anche il nome di “sincronizzazione parziale”

attribuito a questo fenomeno.

Infine, dalle analisi numeriche svolte sul modello originale soggetto a rumore

additivo, uniformemente distribuito e a media nulla, vediamo (figura 2.13(a)) che il

comportamento collettivo persiste, nonostante sia depresso dall’azione del rumore.

Il campo E ha un’ampiezza piu piccola in presenza del rumore, mentre la larghezza

della curva sembra descrescere all’aumentare del numero di oscillatori presenti nella

rete con un andamento del tipo 1/√

N .

Figura 2.13. (a) Rappresentazione nel piano delle fasi (E,E); i punti corrispon-dono a simulazioni effettuate in presenza di rumore additivo uniformemente dis-tribuito. (b) Evoluzione, in funzione del tempo, dell’intervallo temporale tra unosparo e l’altro del singolo neurone (T ). In figura e rappresentato anche un metodoper ricostruire il tempo di sparo del singolo neurone. (c) Evoluzione temporale del

flusso istantaneo π.[22]

2.5 Desincronizzazione in reti neurali diluite

In questa sezione vogliamo analizzare l’effetto della diluizione delle connessioni fra i

neuroni sulla dinamica di una rete neurale con accoppiamenti solo inibitori. Per rete

diluita si intende una rete neurale in cui siano stati recisi dei legami fra neuroni anche

in modo asimmetrico. Per comprendere il ruolo della diluizione, ricordiamo innanzi-

tutto i risultati noti per reti di neuroni LIF globalmente accoppiati. In quest’ambito

sono di notevole importanza i risultati rigorosi pubblicati da Jin in Ref. [30] per una

rete non omogenea di neuroni LIF globalmente accoppiati, le cui connessioni sono

354 356 358 360 362 364Tempo

0

0.5

1

1.5

2

2.5

3

E

Figura 2.14. Andamento del campo E in funzione del tempo nel regime di sin-cronizzazione parziale e in assenza di rumore additivo.

dominate da inibizione globale. In questo caso sia il numero dei neuroni, sia la dis-

tribuzione delle correnti esterne che le ampiezze delle connettivita sono arbitrarie;

ogni neurone inoltre puo avere una diversa soglia di sparo, un diverso potenziale

di riposo della membrana e un diverso valore a cui riassegnare il potenziale dopo

l’emissione di un impulso. Ciononostante Jin ha dimostrato che, in questo sistema,

la sequenza di impulsi converge a configurazioni periodiche stabili da quasi tutti gli

stati iniziali. Ha dimostrato inoltre che il tempo necessario a tale convergenza (detto

transiente) e finito ed al massimo diverge come una potenza del numero di neuroni

nella rete; il transiente si abbrevia al crescere dell’inibizione globale. Si potrebbe

osservare una divergenza del transiente solo nell’eventualita che due neuroni rag-

giungessero contemporaneamente la soglia, ma tale situazione si ottiene solo in un

insieme ristretto di stati iniziali dei neuroni che si restringe sempre di piu al crescere

del livello di inibizione.

A partire dal modello di rete neurale con accoppiamento inibitorio, descritto in

[30], introduciamo la diluizione e mostriamo che nella versione diluita di tale modello

sono osservati tipicamente dei transienti esponenzialmente lunghi, in presenza di un

accoppiamento sufficientemente grande, anche se l’esponente di Lyapunov1 massimo

1La caratterizzazione della dinamica in termini di esponenti di Lyapunov verraa discussa in

del sistema rimane negativo. In altri termini, per piccoli valori della costante di

accoppiamento la dinamica converge, dopo un breve transiente, verso lo splay state

e l’unico effetto dovuto alla soppressione di una piccola percentuale di legami tra

i neuroni e la presenza di inomogeneita nelle mutue interazioni che porta a piccole

difformita negli intervalli tra un’emissione e l’altra dell’impulso. D’altra parte, per

valori dell’accoppiamento sufficientemente grandi, sono osservati dei transienti di

tipo stocastico, la cui durata e esponenzialmente lunga con la dimensione della

rete. Dal momento che vari indicatori mostrano che questo regime e stazionario,

e logico aspettarsi che, nel limite di reti infinitamente grandi, questo transiente

rappresenti una fase termodinamica perfettamente legittima. Dal momento che il

comportamento pseudocaotico e sostenuto dalle discontinuita nella regola con cui

viene effettuata una mappatura [27], e naturale aspettarsi che questo valga anche

nel caso presente; osserviamo infatti che in presenza di disordine, dove i neuroni

non sono piu l’uno equivalente all’altro, si verificano cambiamenti nell’ordine di

sparo che accompagnano le discontinuita nella regola di evoluzione della mappa. La

transizione si manifesta dunque come un fenomeno di desincronizzazione collettivo;

tale fenomeno e analizzabile solo attraverso simulazioni numeriche, dal momento che

l’analisi di stabilita lineare non permette di identificare la soglia oltre la quale si ha

la transizione poiche tutte le traiettorie sono asintoticamente stabili in entrambi i

regimi. Riportiamo di seguito i principali risultati conseguiti da Zillmer et al. [25]

grazie all’analisi numerica.

In una rete di N neuroni LIF, lo stato dell’i-esimo neurone e completamente

determinato dal potenziale di membrana Vi(t) ed obbedisce all’equazione

τ Vi = C − Vi − τ(Vi + W )N∑

j=1

∑

m

gij δ(t − t(m)j ) , (2.52)

dove τ e la costante di tempo di membrana, C e la corrente in ingresso sopra soglia

dettaglio nel capitolo successivo.

(che si riferisce ad una resistenza di membrana unitaria) e W e il potenziale di

inversione. Ogniqualvolta il potenziale Vj(t) raggiunge il valore di soglia Θ, viene

riassegnato al valore R < Θ e viene emesso un impulso; tale impulso e ricevuto

istantaneamente da tutti i neuroni connessi a quello che ha sparato, al tempo t(m)j

(m identifica l’evento di sparo del j-esimo neurone). La ricezione dell’impulso da

parte dell’i-esimo neurone provoca un abbassamento del potenziale di membrana di

tale neurone, in accordo con la trasformazione

V ′i + W = (Vi + W ) exp(−gij). (2.53)

L’ultimo ingrediente che definisce la dinamica del sistema e la matrice di connettivita

gij. La forza di accoppiamento e riscalata con la connettivita del neurone che riceve

l’impulso, per cui

gij =

G/ℓi, se i e j sono connessi,

0, altrimenti,

(2.54)

dove G e la costante di accoppiamento e li il numero di connessioni in ingresso al

neurone i. In altre parole, consideriamo il tipo piu semplice di disordine, determi-

nato dalla presenza, o dall’assenza, dei legami tra i neuroni (le autointerazioni sono

escluse). La frazione rm di legami mancanti e determinata fissando deterministica-

mente il numero totale Nm di legami recisi [Nm = rmN(N−1)]; scegliendo di tagliare

il 5% dei legami si ha rm = 0.05 (caso di debole disordine). A differenza del modello

analizzato in Ref. [30], dove sono distribuite in maniera casuale sia le costanti di

accoppiamento gij, che le correnti di soglia Cij, qui l’unica sorgente di disordine e la

presenza (assenza) di connessioni inibitorie. L’altra importante differenza riguarda

la costante di accoppiamento: mentre in [30] non e ipotizzata nessuna dipendenza

dalla dimensione del sistema, in questo caso e assunta una dipendenza inversamente

proporzionale dal numero di connessioni entranti.

0 0.05 0.1 0.15 0.2r

c

-4.4

-4.2

-4

-3.8

ΛN=100N=200N=300

Figura 2.15. Esponente di Lyapunov massimo in funzione della frazione di legamirecisi e per tre valori diversi della dimensione della rete.[25]

0.5 1 1.5g

103

104

105

<T

r>i,d

N=100N=200

40 60 80 100 120 140N

102

103

104

105

106

107

<T

r>i,d

, <

Tp>

i,d

TransientePeriodoTransientePeriodo

g=1.8

g=2.5

Figura 2.16. (a) Lunghezza media del transiente < ttr > in funzione della costantedi accoppiamento G per due dimensioni diverse della rete. (b) Lunghezza media deltransiente < ttr > (linea continua) e corrispondente periodo < T > dell’attrattoreperiodico (linea tratteggiata) in funzione della dimensione della rete. In ascissa la

scala usata e logaritmica. [25]

Ci aspettiamo che il debole disordine introdotto dalla diluizione riduca la forza

della costante di accoppiamento tra i neuroni e faccia aumentare il valore dell’espo-

nente di Lyapunov massimo Λ, rispetto al caso di neuroni completamente accop-

piati. Analizzando la figura (2.15) vediamo che l’esponente di Lyapunov massimo

cresce al crescere della frazione di legami tagliati, ma rimane negativo almeno fino

a rm = 0.2. Inoltre le simulazioni numeriche indicano che Λ rimane finito nel limite

termodinamico N → ∞.

Dunque, prima o poi, la dinamica deve convergere verso un’orbita periodica

stabile.

La durata del transiente ttr e definita come il tempo piu piccolo necessario affinche

la presente configurazione del potenziale di membrana sia ε vicina ad uno stato

precedente,

ttr = min

t∣

∣

∣

∣

maxi

|vi(t + T ) − vi(t)| < ε

,

dove vi = (Vi − R)/(Θ − R), 1 ≤ i ≤ N e 1 ≤ T ≤ t − 1. Con questa procedura

determiniamo non solo la lunghezza del transiente, ma anche la periodicita T della

soluzione asintotica. Calcoliamo a questo punto la lunghezza media del transiente

< ttr >: la media e effettuata sia su differenti realizzazioni della rete sia su dif-

ferenti condizioni iniziali dei potenziali di membrana. In figura (2.16a) e riportato

l’andamento della lunghezza media del transiente in funzione di G, per due reti di

dimensione diversa (N = 100,200). Per G ≤ 1, < ttr > decresce leggermente al-

l’aumentare della forza di accoppiamento; aumentando la forza di accoppiamento si

verifica un’improvvisa crescita di < ttr >. Inoltre il tasso di crescita e significativa-

mente piu rilevante quando la dimensione della rete viene raddoppiata. Questi dati

suggeriscono dunque l’esistenza di due fasi distinte, che corrispondono approssima-

tivamente a G minore o maggiore di 1. Vediamo adesso come la lunghezza media del

transiente < ttr > e il periodo medio dello stato finale < T > scalano con N (figura

2.16b): per G = 0.5 la lunghezza del transiente cresce linearmente con N , mentre

il periodo medio dell’attrattore asintotico rimane pressoche costante; per G = 1.8

e per G = 2.5 il transiente cresce in modo esponenzialmente veloce. Inoltre per

G = 2.5 anche il periodo cresce esponenzialmente e cio puo essere considerato come

un’indicazione preliminare del fatto che ci possa essere un’ulteriore transizione per

valori della forza di accoppiamento maggiori. Ad ogni modo, focalizziamo la nos-

tra attenzione sui cambiamenti qualitativi che avvengono a G ≈ 1 dove comincia a

crescere esponenzialmente solo il transiente.

10250 10300 10350 10400 10450Tempo

0

10

20

30

40

50

Indi

ce d

i Neu

rone

Figura 2.17. Configurazione di sparo (indice del neurone che spara vs tempo) diun tipico attrattore periodico per G = 2, N = 50 nella LP.[25]

140 150 160 170Tempo

0

50

100

150

200

250

Indi

ce d

i Neu

rone

17350 17400 17450 17500 17550 17600Tempo

0

10

20

30

40

50

Indi

ce d

i Neu

rone

Figura 2.18. (a) Configurazione di sparo per G = 2, N = 50 nella UP. (b) Con-figurazione di sparo associata ad un attrattore periodico nella UP per G = 20,

N = 50. Il periodo e indicato dalle linee verticali.[25]

Per G <≃ 1, l’evoluzione dinamica del sistema e simile a quella del caso com-

pletamente accoppiato. Dopo un transiente di breve durata, il sistema converge ad

uno stato caratterizzato da una sequenza di N spari ( emessi da N neuroni ), che

si ripetono periodicamente (vedi figura 2.17). Tutti i neuroni sparano con lo stesso

ritmo, ma le loro fasi sono equispaziate; in altri termini la soluzione asintotica e

uno stato sincronizzato (in inglese locked phase o LP). La differenza sostanziale con

il caso completamente accoppiato sta nel fatto che le sequenze di impulsi differenti

non sono equivalenti l’una all’altra. Nel limite di neuroni identicamente accoppiati

la dinamica e invariante sotto qualsiasi permutazione, cosa che non e piu vera nel

momento in cui si introducono delle eterogeneita nella connettivita della rete. Ne

deriva cosı che, per rm >≃ 0 esiste un numero esponenzialmente grande di at-

trattori periodici. Piu precisamente, il numero totale di attrattori e dell’ordine di

N !/(∏M

d=1 Nd!), dove Nd e il numero di neuroni connessi agli stessi neuroni (sia in

entrata che in uscita) nella rete ed M e il numero di queste classi di equivalenza

presenti in una data realizzazione della rete.

Per valori della costante di accoppiamento piu grandi, si manifesta un regime

dinamico irregolare (detto unlocked phase o UP), dove il mutuo ordinamento cambia

durante l’intero transiente, come si puo vedere in figura 2.18, dove compare un lento

ma sistematico aggiustamento delle configurazioni.

Riassumendo, la dinamica di transiente ossevata per G ≈> 1: (i) e caratterizzata

da un esponente di Lyapunov negativo; (ii) e effettivamente stazionaria; (iii) ha

durate temporali esponenzialmente lunghe. Queste sono le caratteristiche distintive

del ”caos stabile”, associato alla presenza di discontinuita nello spazio delle fasi.

Nel presente caso la discontinuita e dovuta alla presenza di legami unidirezionali tra

neuroni: nel caso in cui uno dei due legami sia tagliato, uno solo dei due inibisce

l’altro e non si ha piu una mutua inibizione. La presenza di discontinuita nello

spazio delle fasi e una condizione per l’instaurarsi di dinamica caotica stabile, ma

il fatto che non si manifestino transienti anomali per piccoli valori della forza di

accoppiamento indica che la condizione e necessaria ma non sufficiente.

Bibliografia

[1] C. Koch, Biophysics of computation, Oxford University Press, New York

(1999).

[2] H.C Tuckwell, Introduction to theoretical neurobiology - Vol. n. 2, Cam-

bridge University Press, New York (1988).

[3] F. Rieke, D. Warland, Rob de Ruyter van Steveninck e W. Bialek, Spikes: ex-

ploring the neural code, Massachusetts Institute of Technology: Cambridge,

Massachusetts (1996).

[4] A.L. Hodgkin, A.F. Huxley e B. Katz, Measurement of current-voltage

relations in the membrane of the giant axon of Loligo, J. Physiol. 116,

424-448 (1952).

[5] A.L. Hodgkin e A.F. Huxley, Currents carried by sodium and potassium

ion throught the membrane of the giant axon of Loligo, J. Physiol. 116,

449-472 (1952).

[6] A.L. Hodgkin e A.F. Huxley, A quantitative description of membrane

current and its application to conduction and excitation in nerve, J.

Physiol. 117, 500-544 (1952).

[7] S.H. Strogatz, Exploring Complex Networks, Nature 410, 268-276 (2001).

[8] J.D. Murray, Mathematical Biology, Springer-Verlag, Berlino (1989).

[9] W. Gerstner e W. Kistler, Spiking Neuron Models, Cambridge University

Press, Cambridge (2002).

[10] S. Ramon y Cayal, Histologie du Systeme Nerveux de l’Homme et des

Vertebre, A. Maloine, Parigi (1909).

[11] J.G. Nicholss, R.A. Martin e B.G. Wallace, Dai Neuroni al cervello,

Zanichelli, Bologna (1997).

58

[12] P. Dayan e L.F. Abbott, Theoretical Neuroscience: Computational and

Mathematical Modeling of Neural Systems, The MIT Press, Cambridge

(2001).

[13] M.N. Shadlen e W.T. Newsome, The variable discharge of cortical neu-

rons: implications for connectivity, computation, and information

coding, J. Neurosci. 18, 3870-3896 (1998).

[14] D.K.S. Smetters, Electrotonic structure and synaptic integration in cor-

tical neurons, Massachusetts Institute of Technology: Cambridge, Massachu-

setts(1995).

[15] S.G. Cull-Candy, R. Miledi e I. Parker, Single glutamate-activated chan-

nels recorded from locust muscle fibres with perfused patch-clamp

electrodes, J. Physiol. 321, 195-210 (1980).

[16] E.R. Kandel, J.H.Schwartz e T.M. Jessell, Principles of neural science,

McGraw-Hill (2000).

[17] D. Johnston e S. Miao-Sin Wu, Foundations of Cellular Neurophysiology,

Massachusetts Institute of Technology: Cambridge, Massachusetts (1995).

[18] K. Huang, Statistical mechanics, John Wiley & Sons, New York (1987).

[19] D.E. Goldman, Potential, impedance and rectification in membranes,

J. Gen. Physiol. 342, 37-60 (1943).

[20] A.L. Hodgkin e B. Katz, The efffect of sodium ions on the electrical

activity of the giant axon of the squid, J. Physiol. 108, 37-77 (1949).

[21] L.J. Mullins e K. Noda, The influence of sodium-free solutions on mem-

brane potential of frog muscle fibers, J. Gen. Physiol. 47, 117-132 (1963).

[22] P.K. Mohanty e A. Politi, A new approach to partial synchronization in

globally coupled rotators, J. Phys. A: Math. Gen. 39, L415-L421 (2006).

[23] A. Beuter, L. Glass, M.C. Mackey, M.S. Titcombe, Nonlinear dynamics in

physiology and medicine, Springer-Verlag, New York (2003).

[24] S.G Lee, A. Neiman e S. Kim, Coherence resonance in a Hodgkin-Huxley

neuron, Phys. Rev. E 57, 3292-3297 (1998).

[25] R. Zillmer, R. Livi, A. Politi e A. Torcini, Desynchronization in diluted

neural networks, Phys. Rev. E 74, 036203 (2006).

[26] R. Zillmer, R. Livi, A. Politi e A. Torcini, Stability of splay state in pulse-

coupled networks, Phys. Rev. E 76, 046102 (2007).

[27] F. Ginelli, R. Livi, A. Politi e A. Torcini, Relationship between directed

percolation and the synchronization transition in spatially extended

systems, Phys. Rev. E 67, 046217 (2003).

[28] L.F. Abbott e C. van Vreeswijk, Asynchronous states in networks of

pulse-coupled oscillators, Phys. Rev. E 48, 1483 (1993).

[29] Referenza da trovare

[30] D.Z. Jin, Fast convergence of spike sequences to periodic patterns in

recurrents networks, Phys. Rev. Lett. 89, 208102 (2002).

[31] Y. Yu, W. Wang, J. Wang e F. Liu, Resonance-enhanced signal detection

and transduction in the Hodgkin-Huxley neuronal systems, Phys. Rev.

E 63, 21907 (2001).

[32] R.B. Stein, The frequency of nerve action potentials gerated by applied

currents, Proc. Roy. Soc. Lond. B 167, 64-86 (1967).

[33] J. Guckenheimer e R.A. Oliva, Chaos in the Hodgkin-Huxley model, SIAM

J. Applied Dynamical Systems Vol.1, N.1, 105-114 (2002).

[34] J. Rinzel, On repetitive activity in nerve, Fed. Proc. 37, 2793-802 (1978).

[35] J. Rinzel e R. Miller, Numerical calculations of stable and unstable peri-

odic solutions to the Hodgkin-Huxley equations, Math. Biosci. 49, 27-59

(1980).

[36] T.B. Kepler, L.F. Abbott e E. Marder, Reduction of conductance-based

neuron models, Biol. Cybern. 66, 381-387 (1992).

[37] B. van der Pol, On relaxation oscillations, Phil. Mag. 2, 978-992 (1926).

[38] J. Cronin, Mathematical Aspects of Hodgkin-Huxley Neural Theory,

Cambridge University Press: Cambridge, Regno Unito (1987).

[39] C.E. Shannon, A mathematical theory of communication, The Bell Sys-

tem Technical J. 27, 623 (1948); 27, 379 (1948).

[40] G. Boffetta, M. Cencini, M. Falcioni e A. Vulpiani, Predictability: a way to

characterize complexity, Physics Reports 356, 367-474 (2002).

[41] D.M. MacKay e W.S. McCulloch, The limiting information capacity of a

neuronal link, Bull. math. Biophysics 14, 127 (1952)

[42] J.D. Farmer, Information Dimension and the Probabilistic Structure

of Chaos, Z. Naturforsch. 37 a, 1304-1325 (1982).

[43] F. Moss, A. Bulsara e M.F. Shlesinger, Proceedings of the NATO Ad-

vanced Research Workshop on Stochastic Resonance in Physics and

Biology, J. Stat. Phys. 70 (1993).

[44] R. Mannella e P.V.E. McClintock, Proceedings of the International Work-

shop on Fluctuations in Physics and Biology: Stochastic Resonance,

Signal Processing and Related Phenomena, Nuovo Cimento Soc. Ital.

Fis. 17 (1995).

[45] L. Gammaitoni, P. Hanggi, P. Jung e F. Marchesoni, Stochastic resonance,

Reviews of Modern Physics 70, 223-287 (1998)

[46] H.E. Plesser e W. Gerstner, Noise in integrate and fire neurons: from

stochastic input to escape rates, Neural Comput. 12, 367-384 (2000).

[47] A. Neiman, B. Shulgin, V. Anishchenko, W. Ebeling, L. Schimansky-Geier e J.

Freund, Dynamical Entropies applied to Stochastic Resonance, Phys.

Rev. Lett. 76, 4299-4302 (1996).

[48] E. Salinas e T.J. Sejnowski, Impact of correlated synaptic input on out-

put firing rate and variability in simple neuronal models, J. Neurosci.

20, 6193-6209 (2000).

[49] D.R. Cox, Renewal Theory, Methuen, Londra (1962).

[50] D.R. Cox e H.D. Miller, The Theory of Stochastic Processes, Wiley, New

York (1965).

[51] D. Brown, J. Feng e S. Feerick, Variability of firing of Hodgkin-Huxley

and FitzHugh-Nagumo neurons with stochastic synaptic input, Phys.

Rev. Lett. 82, 4731-4734 (1999).

[52] J. Feng, Behavior of integrate-and-fire and Hodgkin-Huxley models

with correlated inputs, Phys. Rev. E. 63, 51902 (2001).

[53] A.S. Pikovsky e J. Kurths, Coherence Resonance in a Noise-Driven Ex-

citable System, Phys. Rev. Lett. 78, 775-778 (1997).

[54] N.G. van Kampen, Stochastic processes in physics and chemistry, North-

Holland personal library, Amsterdam (1997).

[55] H.A. Kramers, Brownian motion in a field of force and the diffusion

model of chemical reaction, Physica 7, 284 (1940).

[56] B. Lindner, J. Garcıa Ojalvo, A. Neiman e L. Schimansky-Geier, Effects of

noise in excitable systems, Physics Reports 392, 321-424 (2004).

[57] W.J. Wilbur e J. Rinzel, A Theoretical Basis for Large Coefficient of

Variation and Bimodality in Neuronal Interspike Interval Distribu-

tions, J. theor. Biol. 105, 345-368 (1983).

[58] L. Sacerdote e R. Sirovich, Multimodality of the interspike interval distri-

bution in a simple jump-diffusion model, Scientiae Mathematicae Japon-

icae Online 8, 359-374 (2003).

[59] G. Giacomelli, M. Giudici, S. Balle e J.R. Tredicce, Experimental evidence

of coherence resonance in an optical system, Phys. Rev. Lett. 84, 3298

(2000).

[60] S. Zhong e H. Xin, Noise-induced oscillations and internal stochastic

resonance in a model of excitable biomembrane, Chem. Phys. Lett. 321,

309 (2000).

[61] A. Ganopolski e S. Rahmstorf, Abrupt glacial climate changes due to

stochastic resonance, Phys. Rev. Lett. 88, 38501 (2002).

[62] A. M. Lacasta, F. Sagues e J.M. Sancho, Coherence and anticoherence

tuned by noise, Phys. Rev. E 66, 45105 (2002).

[63] E. Salinas e T.J. Sejnowski, Correlated neuronal activity and the flow of

neuronal information, Nature Reviews in Neuroscience 2, 539 (2001).

[64] P. Berge, Y. Pomeau e C. Vidal, Order within chaos, Herman e John Wiley

& Sons, Parigi (1984).

[65] J. Clinton Sprott, Chaos and Time-Series Analysis, Oxford University

Press, New York (2003).

[66] W.H. Press, S.A. Teukolsky, W.T. Vetterling e B. P. Flannery, Numerical

recipes in fortran, Cambridge University Press (1992).

[67] Th. Kreuz, comunicazione privata (2005).

Dinamica di reti neuronali - Istituto Nazionale di Fisica ... 2.4.1 Lunghezze d’onda grandi ... un...

Documents

Transcript of Dinamica di reti neuronali - Istituto Nazionale di Fisica ... 2.4.1 Lunghezze d’onda grandi ... un...