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