Figures for Gonzalo’s talk in Tarragona

>> load(‘c:\hipertec\AblationElegans\Datos\datos_celegans_paradinamica.mat’)
>> V0=zeros(279,1);
V0([40 41])=1;
>> M=conectividad2matrizsistema(todas.A_chem,todas.A_ej,todas.GABA,1,1,30);
>> V=matrizsistema2V_num(M,V0,[-1 1],.0005,1000);
>> imagesc(V)
>> caxis([0 .3])

>> resultados=matrizsistema2experimentoablacion(M,V0,[-1 1],54:57,.0005,1000);
>> load(‘c:\hipertec\AblationElegans\Datos\info_redes.mat’)
>> [roc,umbral]=resultados2pseudoroc(resultados,chemotaxis1.buenas,[chemotaxis1.sensiales chemotaxis1.finales],1);

Good.

>> find(roc(:,2)>.8,1)
ans =
36
>> umbral(36)
ans =
0.0017

>> redencontrada=find(mean(resultados,2)>.0017);
>> pintared_posreal(todas.A_chem,todas.A_ej,todas.nombres,1:279,chemotaxis1.sensiales,chemotaxis1.finales,redencontrada)

Advertisement

Saving information about known networks

>> chemotaxis1.info=’Red de chemotaxis que aparece en el primer post sobre dinámica de elegans’;
>> chemotaxis1.buenas=buenas;
>> chemotaxis1.buenas=[11:14 19:20 40:41 54:57 77:78 111:113 168:169 178:179];
>> chemotaxis1.sensoriales=[40 41];
>> chemotaxis1.motoras=54:57;
>> save info_redes chemotaxis1

 

Comparison of different dynamic models for C elegans

I use something like a ROC curve, that measures the proportion of the neurons belonging to the subnetwork and the proportion of neurons not belonging to the subnetwork that the method finds, depending on the threshold:

I use the chemotaxis network I used before.

Comparison between the model with all conductances equal (as in the older post) and the model with all connections as if gap junctions, and no membrane conductance:

>> clear
>> load(‘c:\hipertec\ablationelegans\Datos\datos_celegans_paradinamica.mat’)
>> V0=zeros(279,1);
>> V0([40 41])=1;
>> M{1}=conectividad2matrizsistema(todas.A_chem,todas.A_ej,todas.GABA,1,1,1);
>> resultados{1}=matrizsistema2experimentoablacion(M{1},V0,500,.005,[-1 1],54:57);
>> M{2}=conectividad2matrizsistema(zeros(279),todas.A_chem+todas.A_ej,todas.GABA,1,1,0);
>> resultados{2}=matrizsistema2experimentoablacion(M{2},V0,500,.005,[-1 1],54:57);
>> buenas=[11:14 19 20 77 78 111:113 168 169 178 179]; % Neurons belonging to the subnetwork
>> resultados2pseudoroc(resultados{1},buenas,[40 41 54:57],’b’);
>> hold on
>> resultados2pseudoroc(resultados{2},buenas,[40 41 54:57],’r’);

 

Now I compare the old model (with all conductances equal) with different time scales:

>> resultados{1,2}=matrizsistema2experimentoablacion(M{1},V0,250,.005,[-1 1],54:57);
>> resultados{1,3}=matrizsistema2experimentoablacion(M{1},V0,100,.005,[-1 1],54:57);
>> clf
>> resultados2pseudoroc(resultados{1,1},buenas,[40 41 54:57],’b’);
>> hold on
>> resultados2pseudoroc(resultados{1,2},buenas,[40 41 54:57],’r’);
>> resultados2pseudoroc(resultados{1,3},buenas,[40 41 54:57],’k’);

Perfect overlapping, so the time scale does not make any effect in this range.

Test of the pseudoROC:

>> resultados=rand(279,1);
>> clf
>> resultados2pseudoroc(resultados,buenas,[40 41 54:57],’k’);

OK.

Now, I compare the model with equal conductances with a model tuned so that it does not saturate:

>> M{1}=conectividad2matrizsistema(todas.A_chem,todas.A_ej,todas.GABA,1,1,1);
>> resultados{1}=matrizsistema2experimentoablacion(M{1},V0,500,.005,[-1 1],54:57);
>> resultados{1}=matrizsistema2experimentoablacion(M{1},V0,500,.005,[-1 1],54:57);
>> M{2}=conectividad2matrizsistema(todas.A_chem,todas.A_ej,todas.GABA,.05,1,1.5);
>> V=matrizsistema2V_num(M{2},V0,10000,.005,[-1 1]);
>> plot(V(54:57,:)’)

>> resultados{2}=matrizsistema2experimentoablacion(M{2},V0,10000,.005,[-1 1],54:57);
>> resultados2pseudoroc(resultados{2},buenas,[40 41 54:57],’r’);
>> hold on
>> resultados2pseudoroc(resultados{1},buenas,[40 41 54:57],’b’);

The model with equal conductances works better.

First results for network identification

I just run the network with the linear model with saturation, and all conductances equal to 1. Then, I ablate neurons one by one, and see the change in the target neurons. I am looking for the following network:

>> M=conectividad2matrizsistema(todas.A_chem,todas.A_ej,todas.GABA,1,1,1);

Initial conditions, the sensory neurons activated:

>> find(strcmpi(todas.nombres,’aser’))
ans =
41
>> find(strcmpi(todas.nombres,’asel’))
ans =
40
>> V0=zeros(279,1);
>> V0([40 41])=1;

Target neurons:

>> find(strcmpi(todas.nombres,’avar’))
ans =
55
>> find(strcmpi(todas.nombres,’aval’))
ans =
54
>> find(strcmpi(todas.nombres,’avbl’))
ans =
56
>> find(strcmpi(todas.nombres,’avbr’))
ans =
57
>> resultados=matrizsistema2experimentoablacion(M,V0,500,.005,[-1 1],54:57);
>> plot(resultados)

The highest peaks are for the target neurons. Zoom:

The high peaks around 40 are the target neurons, but the other high peak are other neurons:

>> todas.nombres(find(mean(resultados,2)>4*10^-3 & mean(resultados,2)<10^-2))
ans =
‘AIBL’
‘AIBR’
‘ASEL’
‘ASER’

AIB, that plays a central role in the network. As we lower the threshold, we find more and more neurons in the network, and also some which are not:

>> todas.nombres(find(mean(resultados,2)>1*10^-3 & mean(resultados,2)<10^-2))
ans =
‘AIBL’
‘AIBR’
‘AIYL’
‘AIYR’
‘AIZL’
‘AIZR’
‘ASEL’
‘ASER’
‘AVEL’
‘AVER’
‘PVCL’
‘PVCR’
‘RIBL’
‘RIML’
‘RIMR’

We catch 7 out of 15 neurons in the network (counting laterally doubled neurons, excluding sensory and target neurons). 6 neurons are not in the network.

>> todas.nombres(find(mean(resultados,2)>.5*10^-3 & mean(resultados,2)<10^-2))
ans =
‘AIAL’
‘AIAR’
‘AIBL’
‘AIBR’
‘AIYL’
‘AIYR’
‘AIZL’
‘AIZR’
‘ASEL’
‘ASER’
‘AVDR’
‘AVEL’
‘AVER’
‘DVA’
‘FLPL’
‘PVCL’
‘PVCR’
‘RIBL’
‘RIBR’
‘RIML’
‘RIMR’
‘SAAVL’
‘VA08’

We catch 11 out of 15 neurons. 10 neurons are not in the network.

>> todas.nombres(find(mean(resultados,2)>.4*10^-3 & mean(resultados,2)<10^-2))
ans =
‘AIAL’
‘AIAR’
‘AIBL’
‘AIBR’
‘AIYL’
‘AIYR’
‘AIZL’
‘AIZR’
‘AS07’
‘ASEL’
‘ASER’
‘AVDL’
‘AVDR’
‘AVEL’
‘AVER’
‘AWCL’
‘DVA’
‘DVC’
‘FLPL’
‘PVCL’
‘PVCR’
‘RIBL’
‘RIBR’
‘RIFR’
‘RIML’
‘RIMR’
‘SAAVL’
‘VA08’
‘VB09’

We catch 13 out of 15 neurons. 14 neurons are not in the network.

 

Excitation of neurons connected to “the nose”

sensoriales_nariz=find(todas.S(:,1));
>> V0=zeros(279,1);
>> V0(sensoriales_nariz)=1;
>> solucion=matrizsistema2solucion(M,V0);
>> V=solucion2voltajes(solucion,0:.001:.5);
>> [s,orden]=sort(todas.pos_real);
>> imagesc(V(orden,:))

Note that neurons do not have the usual order. They are ordered by position in the worm.

>> imagesc(log10(abs(V(orden,:))))

Longer times:

>> V=solucion2voltajes(solucion,0:.005:2);

I do not see any propagation…

First results for the dynamics of C elegans

Eigenvectors and eigenvalues MUST be complex. Still, real solution for the voltage is possible (each eigenvector/value has another one which is complex conjugate. If the modulus of the coefficient of both is the same, the final linear combination is real).

With no membrane conductance:

>> M=conectividad2matrizsistema(todas.A_chem’,todas.A_ej,todas.GABA,1,1,0);
>> V0=zeros(279,1);
>> solucion=matrizsistema2solucion(M,V0);
>> plot(solucion.autovalores,’.’)

This figure matches with Figure 8 of Mitya’s paper (y axis is imaginary part)

We see that there are positive eigenvalues. This has two consequences:
– Diverging voltage.
– Apparently some computations are not performed correctly. I get non-negligible imaginary parts in the voltages.

I see that I need about 30 in the membrane conductance to get all eigenvalues negative:

>> M=conectividad2matrizsistema(todas.A_chem’,todas.A_ej,todas.GABA,1,1,30);
>> solucion=matrizsistema2solucion(M,V0);
>> plot(solucion.autovalores,’.’)

 

FIRST EXAMPLE: Excitation of neuron 1.

>> V0=zeros(279,1);
>> V0(1)=1;
>> solucion=matrizsistema2solucion(M,V0);
>> V=solucion2voltajes(solucion,0:.001:.5);
>> plot(V(1,:))

>> plot(V(2,:))

And so on…

SECOND EXAMPLE: Random excitation

>> V0=rand(279,1);
>> solucion=matrizsistema2solucion(M,V0);
>> V=solucion2voltajes(solucion,0:.001:.5);
>> imagesc(V)

What are those ‘rough’ zones?

>> plot(V(149,:))

What the…?

Addition of GABA-ergic neurons to the C elegans data

I use the data given by Varshney et al. in http://mit.edu/lrv/www/elegans/

>> load(‘C:\Hipertec\ablationelegans\datos\Datos_VarshneyEtAl\ConnOrdered_040903(2).mat’)
for c=1:279
ind(c)=find(strcmpi(todas.nombres{c},Neuron_ordered));
end    
>> todas.GABA=logical(GABAergic(ind));
>> todas.nombres(todas.GABA)
ans =
‘AVL’
‘DD01’
‘DD02’
‘DD03’
‘DD04’
‘DD05’
‘DD06’
‘DVB’
‘RIS’
‘RMED’
‘RMEL’
‘RMER’
‘RMEV’
‘VD01’
‘VD02’
‘VD03’
‘VD04’
‘VD05’
‘VD06’
‘VD07’
‘VD08’
‘VD09’
‘VD10’
‘VD11’
‘VD12’
‘VD13’

I check in the original reference, McIntire SL, Jorgensen E, Kaplan J, Horvitz HR (1993) The GABAergic nervous system of Caenorhabditis elegans. Nature 364:337–341. doi:10.1038/364337a0, that these are indeed the 26 identified neurons that belong to the GABAergic system in C elegans.

 

>> save datos_celegans_paradinamica todas

New compilation of C elegans data from the original Excel files

The old data (stored in datos_celegans_struct) have a symmetrical connectivity matrix. I repeat the extraction from the original Excel files, with a new program:

>> todas=cogedatos_celegans_paradinamica(‘NeuronConnect.xls’,’NeuronFixedPoints.xls’,’NeuronType.xls’)

Importando datos desde archivos de Excel…Terminado
Nombres y posiciones extraidos: 279 neuronas encontradas

Extrayendo matriz de conectividad neurona-neurona (A)…
1000 entradas leidas
2000 entradas leidas
3000 entradas leidas
4000 entradas leidas
5000 entradas leidas
6000 entradas leidas
Neurona VC06 no encontrada. Conecta con NMJ.
Finalizado: 6417 entradas leidas. Matriz A generada.
Comprobación de simetría de A_ej: Correcto.
Comprobación de sinapsis químicas enviadas = sinapsis químicas recibidas: Correcto.
Conexiones eléctricas: 1777.
Conexiones químicas únicas: 2430.
Conexiones químicas múltiples: 3964.
Uniones neuromusculares: 1409.

Extrayendo matrices de conectividad neurona-organo (S y M)…
Hay más uniones neuromusculares en NeuronConnect que conexiones neurona-músculo en NeuronFixedPoints para las siguientes neuronas:
ans =
‘PLNL’
‘VA09’
‘VB10’

NOTA: El 6 de junio de 2009 comprobé que esto se da para las neuronas PLNL, VA09 y VB10 en
los archivos de Excel originales. Para estas neuronas, no es un error del programa de extacción.
Finalizado: 650 entradas leidas. Matrices S y M generadas, con 108 órganos.
Conexiones a sensores: 86.
Conexiones a músculos: 1865.1.

todas =
nombres: {279×1 cell}
pos_real: [279×1 double]
ganglios: [279×1 double]
nombres_ganglios: {1×10 cell}
ambiguity: {279×1 cell}
cod_amb: [279×1 double]
cod_amb2ambiguity: {‘None’ ‘MB’ ‘MTS’ ‘MAS’ ‘MD’ ‘MAD’ ‘MS’ ‘RDI’ ‘RDM’ ‘RVI’}
A_ej: [279×279 double]
A_chem: [279×279 double]
S: [279×108 double]
M: [279×108 double]
f: [1×108 double]
nombres_organos: {1×108 cell}

 

NOTE: The argument factor_poly was 1 by default. Therefore, no distinction is made between single chemical synapses and multiple chemical synapses.

DIFFERENCES WITH THE OLD DATA:
– Chemical synapses were repeated, to symmetrize.
– Electrical synapses were divided by two, in order to count them only once (maybe not so good idea, because the equation for the cost already takes this into account?). Here, I keep them as they are.