20 mai 2008 – analyse des résultat de sage – restriction WT et double mutant (suite)

10h40
après avoir obtenu les segment, je vais donc chercher a les localiser par rapport aux orf environnantes.
Tout d’abord j’exporte la liste des segment sur douanierousseau
scp berthe:sage/total/stats/ratio2_5/WT&mutant/segments_WT_CUTs ~/sage/total/stats/ratio2_5/WT&mutant/
ensuite j’utilise segment.py (dans le workspace genepy) pour transformer ce fichier du format matlab au format csv puis pour localiser les orf
j’ai modifier légèrement segment.py pour que les nom des colonnes correspondent bien a ce qui est dedans.

j’obtiens le fichier segments_WT_CUTs_density_max.csv que j’exporte sur bambino pour appliquer la macro excel de search_features_A_and_B.xls dessus
j’enregistre les résultats dans segments_WT_CUTs_density_max.xls et j’exporte dans le fichier segments_WT_CUTs_dist_to_feature.txt
j’élimine les colonne de texte et j’enregistre en csv
puis j’exporte sur berthemorisot et je l’importe dans matlab

j’obtiens la matrice seg_WT_CUTs_dist2feature et la matrice colheader_seg_WT_CUTs_dist2feature qui contient les intitulés des colonnes.
je représentes les nb de segment de chaque type en fontion de leur distance respectives avec soit le debut ou la fin des orf
j’ai écrit un script qui traces les courbes de façon automatique en
mettant la plupart des bon parametres

function [distance_hist]=plot_dist2feature(data, n, smooth_window)

extreme=cell(1,2);
extreme{1,1}='start';
extreme{1,2}='stop';

sens=cell(1,2);
sens{1,1}='sens';
sens{1,2}='antisens';

classe=cell(2,2);
classe{1,1}='Oligo-dT';
classe{2,1}=8:10;
classe{1,2}='CUTs';
classe{2,2}=1:2;

cols=cell(4,3);

cols{1,1}='start';
cols{1,2}='sens';
cols{1,3}=[8 10 16];

cols{2,1}='start';
cols{2,2}='antisens';
cols{2,3}=[12 14 18];

cols{3,1}='stop';
cols{3,2}='sens';
cols{3,3}=[9 11 17];

cols{4,1}='stop';
cols{4,2}='antisens';
cols{4,3}=[13 15 19];

max_plot=0;
min_plot=0;
distance_hist=cell(8,5);
i=1;
for id_classe=1:2
for id_extreme=1:2
for id_sens=1:2
for j=1:4
feature_extreme=extreme{1,id_extreme};
feature_sens=sens{1,id_sens};
if strcmp(cols{j,1}, feature_extreme) && strcmp(cols{j,2}, feature_sens)
columns=cols{j,3};
end
end

distance_hist{i,1}=mk_distance_hist(data, classe{2,id_classe}, columns, -1000, 1000, n, 20);
distance_hist{i,2}=feature_extreme;
distance_hist{i,3}=feature_sens;
distance_hist{i,4}=classe{1,id_classe};
distance_hist{i,5}=n;
h=zeros(2001,1);
tmp = distance_hist{i,1};
h(min(tmp)+1001:max(tmp)+1001)=hist(tmp,max(tmp)-min(tmp)+1)';
distance_hist{i,6}=h;
distance_hist{i,7}=smooth(h, smooth_window, 0);
i=i+1;
end
end
end
zero_line=[0 0];
figure
set(gcf,'Color','k');
plot_ids = 1:2:7;
for i=1:4
if i max_plot
max_plot= limite(2);
end
end

for i=1:4
subplot(2,2,i)
ylim([min_plot max_plot]);
line(xlim, zero_line, 'Color', 'w');
line(zero_line, ylim, 'Color', 'w');
end


j’ai rajouté qq éléments sur le graphique.
ce qui me donne ca pour n=4 :

quelques chiffres sur ces données
on a sélectioné uniquement les tag provennant de la souche sauvage (oligo-dT) et du double mutant Trf4/Rrp6 (CUTs)
on a :

oligo-dT -> 49284 tags
Cuts -> 66884 tags

environ 1/3 des tags des Cuts a été éliminé de façon aléatoire afin d’avoir le meme nombre pour les 2 catégories.

on a alors

oligo-dT -> 49284 tags
Cuts -> 49077 tags

il y a au total 98361 tags

la détection des segments sur ces tags
permet de trouver 12021 clusters
ces clusters contiennent au moins 1 tag
si on augmente le nb minimal de tag on a :

1-> 12021
2-> 6001
3-> 4047
4-> 3069
5-> 2482
6-> 2072
7-> 1791
8-> 1572
9-> 1410
10-> 1277

le calcul du ratio Oligo-dT / (CUTs + Oligo-dT) puis le fractionnement des clusters dans 10 classes en fonction de ce ratio donne

Nb tag minimum	Classes -> 1	2	3	4	5	6	7	8	9	10
1	6168	179	94	178	330	60	161	216	200	4435
2	2573	179	94	178	330	60	161	216	200	2010
3	1568	179	94	178	97	60	161	216	200	1294
4	1169	179	94	62	97	60	59	216	200	933
5	935	179	46	62	58	60	59	148	200	735
6	787	140	46	42	58	30	59	107	200	603
7	691	117	46	27	40	30	46	107	176	511
8	621	100	37	27	32	20	46	97	152	440
9	572	85	28	24	29	20	41	85	137	389
10	538	73	19	18	24	16	37	74	121	357

on va considerer que les CUTs sont les segments des classes 1 et 2 et les oligo-dT ceux des classes 8,9 et 10. de plus on ne va s’intéresser qu’aux segments ayant au moins 4 tags
ce qui fait
CUTs -> 1348 segments
oligo-dT -> 1349 segments

dans ces segments on a
38 % (508 segments) des oligo-dT sont entre 600 et 0 avant le start des orf en sens
88 % (1198 segments) des oligo-dT sont entre 0 et 400 après le stop des orf en sens
48 % (646 segments) des oligo-dT sont entre -100 et 400 par rapport au stop des orf en antisens
45 % (605 segments) des CUTs sont entre -600 et 400 par rapport au start des orf en sens
45 % (613 segments) des CUTs sont entre -1000 et -200 par rapport au start des orf en antisens
40 % (538 segments) des CUTs sont entre -600 et 100 par rapport au stop des orf en antisens

15 mai 2008 – analyse des résultat de sage – restriction WT et double mutant

10h25
Je vais reprendre l’analyse des sage a la base en ne prennant en compte que les tag obtenu sur les manip wild type et double mutant Trf4/Rrp6

Je dois refaire le clustering afin de n’avoir que des segment correspondant a des groupes de tag de ces manip.

1- ouverture dans matlab du fichier
~/sage/total/stats/tag_tot_genome.mat

j’ai dans ce fichier les matrice correspondant aux 12 manip de sage (header01 à header12)
chaque matrice à une taille de 24 million ce qui corresponde au nombre de nucleotide de l’ensemble du génome mis bout a bout en séparant les brin C et W
les matrice chr_len et chr_ref donnent les positions de chaque chromosome dans cette matrice

2- la recherche des segments (cluster) va se faire uniquement sur les manip 1 et 8, correspondant au Wild Type (WT) et au double mutant Trf4/Rrp6

3- j’utilise le programme segment.m avec la commande
header = header01 + header08;
header(:,2) = header01;
[header_region] = segment(header, dt, nt, chr_ref);
je récupère la matrice ratio2_1_100_1 qui a été crée
et je représente les histogramme des classe de ratio

figurefor i=1:9subplot(3,3,i)hist(ratio2_1_100_1(find(ratio2_1_100_1(:,4)>i-1),1))xlabel('Ratio Oligo-dT / (CUT + Oligo-dT)', 'color',[1 1 1]);ylabel('Number of clusters', 'color',[1 1 1]);title(strcat('n = ', num2str(i)), 'color',[1 1 1])set(gca, 'Xcolor', [1 1 1], 'Ycolor', [1 1 1])end

et j’obtiens la figure suivante

on voit qu’il y a un nb plus important de cluster de la classe 1 que ceux de la classe 10
si on regarde le nb de tag, on constate aussi qu’il y a environ 20 % de plus de CUTSque d’Oligo-dT

On va donc chercher a réduire artificiellement le nb de CUT en en éliminant de façon aléatoire jusqu’à avoir le meme nombre
script : delete_rand_tag

h8 = header08;tag = find(h8);s = size(tag,1);while sum(h8) > 49284p = floor(rand(1000,1) * s);ptag = tag(p);h8(ptag) = 0;end

header_rand = header01 + h8;
header_rand(:,2) = header01;
je fais ensuite [header_region_rand] = segment(header_rand, dt, nt, chr_ref);
je renome ratio2_1_100_1 en ratio2_1_100_1_ rand
puis je représente les memes histogrammes

figurefor i=1:9subplot(3,3,i)hist(ratio2_1_100_1_rand(find(ratio2_1_100_1_rand(:,4)>i-1),1))xlabel('Ratio Oligo-dT / (CUT + Oligo-dT)', 'color',[1 1 1]);ylabel('Number of clusters', 'color',[1 1 1]);title(strcat('n = ', num2str(i)), 'color',[1 1 1])set(gca, 'Xcolor', [1 1 1], 'Ycolor', [1 1 1])end

et j’obtiens

on voit que c’est plus homogène entre les classes extrèmes.
pour vérifier cela, on cherche a voir si la répartition du nombre de tag dans les segment est la meme pour toute les classes.

j’extrait de la matrice ratio2_1_100_1_ rand les lignes pour lesquelles le ratio est faible (inf a 0.2 CUTs) et fort (sup à 0.8 WT) et je fait l’histogramme du nombre de tag (colonne 4). je choisi pour ces histogramme d’avaoir comme nombre de catégorie la valeur maximale de tag pour chaque type (3318 pour le WT et 1522 pour les CUTs.
puis je trace les courbe correspondantes en échelle log

h_0=hist(ratio2_1_100_1_rand(find(ratio2_1_100_1_rand(:,1)0.8),4), max(ratio2_1_100_1_rand(find(ratio2_1_100_1_rand(:,1)>0.8),4)));
figure;
plot(h_0,'r');
hold on
plot(h_1,'g');

on voit que les segments avec un faible nb de tag sont plus nombreux pour les CUTs. cependant a partir de 4 tags par segment on observe la meme répartition.
on peut donc dire que a partir de 4 tags par segment on a un nombre similaire de segments avec le meme nombre de tag dans ceux-ci.

Je me suis rendu compte que les position des tags sont fausses. en fait elles sont décalée par rapport a ce qu’on veux car on s’interresse a la fin des tag et non au début
j’ai fait un programme qui corrige ca (tag_position_correction)
je corrige donc les 2 matrice contenant les résultats des WT et CUTs

chr_WT=tag_position_correction(header01,chr_len);chr_CUTs=tag_position_correction(h8,chr_len);chr_WT_CUTs = chr_WT + chr_CUTs;chr_WT_CUTs(:,2)=chr_WT;[regionWT_CUTs] = segment(chr_WT_CUTs, dt, nt, chr_ref);

cpt =

       12021

et je récupère la matrice segments_WT_CUTs
je vais maintenant rechercher le max de densité de chaque segment
je calcul d’abord la matrice de densité pour les tag WT et CUts
density_WT_CUTs = densite_sage(chr_WT_CUTs(:,1), 100);
puis je cherche le max de densité
segments_WT_CUTs=find_max_segment(density_WT_CUTs, segments_WT_CUTs);
j’exporte avec save(strcat('segments_WT_CUTs','.mat.csv'),'segments_WT_CUTs', '-ascii','-tabs');

TODO -> export sur douanierousseau, localiser les segment par rapport aux orf et faire les graph