view src/mz_4_extract_final_mz_v3.sci @ 2:05a5f68935f2 draft default tip

Deleted selected files
author jcb-mpl
date Thu, 12 May 2022 12:14:45 +0000
parents b2d949399888
children
line wrap: on
line source

function [list_eic4_files,mz4_stats,mz3]=mz_4_extract_final_mz_v3(list_eic_files,delta_mz)


    // list_eic_files: une liste de fichiers 
    //                 chaque fichier = (TR x mz); variable =  x_tr_mz3 
    // mz_uncertainty: l'incertitude sur les m/z utilisee pour calculer les matrices d'EICs; typiquement: 0,0024 
    // mz3_stats: donne les moyennes min et max du même m/z obtenues sur les differents fichiers 
    //            et calcule l'intervalle pour reconstruire les matrices d'EICs 
    // -------------------------------------------------------------------------
    
    // regrouper tous les m/z dans un seul fichier 
    // colonne 1 = valeur de m/z 
    // colonne 2 = code du fichier; 1 à 6 si 6 fichiers 
    // colonne 3 = n° de ligne dans chaque fichier 
    
    // ex de list_eic_files: 
    //     list_eic_files=[ 'ID1_EICs-raw.mat'; 'ID2_EICs-raw.mat';'ID3_EICs-raw.mat';'AC-ID1_EICs-raw.mat';'AC-ID2_EICs-raw.mat';'AC-ID3_EICs-raw.mat'];
    
    if argn(2)<2 then
        delta_mz=0.0024;
    end

    n_files=max(size(list_eic_files));
    
    mzs=[];   // mzs=compilation de tous les mz
    
    for i=1:n_files;
        x_tr_mz3=list_eic_files(i); 
 
        list_mz=strtod(x_tr_mz3.v);
        n=size(list_mz,1);  
        if i==1 then 
            mzs=[list_mz ones(n,1) [1:n]' ];
        else
            mzs=[mzs; [list_mz i*ones(n,1) [1:n]']]; 
        end
    end    
   
    // tri selon des valeurs croissantes 
    [nul, index]=gsort(mzs(:,1),'g','i');
    mzs=mzs(index,:);
    
    
    n=size(mzs,1); 
    mz2=zeros(n,n_files,2);       // 6 remplace par n_files échantillons   12dec19
                                  // 2 colonnes: mz et n° de pic
    
    // choix de l'incertitude 
    // delta_mz=0.0024; OK 
    // delta_mz=0.0012; essayé mais pas concluent: des m/z se trouvent décalés seuls sur une ligne 
    
    index=1;
    for i=1:n-1;
        if mzs(i+1,1)-mzs(i,1)>delta_mz | (mzs(i+1,1)-mzs(i,1)<= delta_mz & mz2(index,mzs(i+1,2),1)~=0 ) then
            index=index+1;
        end
        mz2(index,mzs(i+1,2),1)=mzs(i+1,1);  
        mz2(index,mzs(i+1,2),2)=mzs(i+1,3);
    end                 
                
    mz2=mz2(1:index,:,:);     // 10683 m/z  
    // visuellement: le fichier semble OK 
    
    // tri des plus intéressants : 3 dans un groupe au moins 

    // commente le 12dec19  => il faut retrouver dans toutes les répétitions
    tri=mz2(:,1,1);
    for k=2:n_files;
        tri=tri.*mz2(:,k,1);
    end
    // tri1=mz2(:,1,1).*mz2(:,2,1).*mz2(:,3,1);
    // tri2=mz2(:,4,1).*mz2(:,5,1).*mz2(:,6,1);
    // tri=tri1+tri2;


    indice=find(tri >0);
    mz3.d=mz2(indice,:,:);       // 8589 m/z
    // visuellement: semble OK
    
    // rajout des statistiques
    n=size(mz3.d,1);
    mz3_min_max=zeros(n,5);
    for i=1:n;
        xtemp=mz3.d(i,1:n_files,1);
        xtemp=matrix(xtemp,[n_files,1,1]);
        index=find(xtemp~=0);
        mztemp=xtemp(index);
        mean_temp=mean(mztemp);

        n_temp=max(size(mztemp));
        std_temp=std(mztemp)*sqrt(n_temp/(n_temp-1));
        mz3_min_max(i,1)=min(xtemp(index));
        mz3_min_max(i,2)=max(xtemp(index));
        mz3_min_max(i,3)=mean_temp;
        mz3_min_max(i,4)=mean_temp - 2*std_temp;
        mz3_min_max(i,5)=mean_temp + 2*std_temp;
    end
    
    // mise en div:
    mz3.i=string(mz3_min_max(:,3));
    mz3.v.v1='file'+string([1:size(mz3.d,2)]');
    mz3.v.v2=['signal';'column'];
//pause
    mz3=div(mz3);
    
    mz3_stats.d=mz3_min_max;
    mz3_stats.i='M'+ string(mz3_min_max(:,3));
    mz3_stats.v=['mean min';'mean max';'mean mz';'mean -2std';'mean +2std'];
    mz3_stats=div(mz3_stats);
   
 
    // --------------------------------------------------------------------------
    // a ce stade on a :
    //           mz3 =  une liste de m/z commune et interessante  (8589 x 6 x 2)
    //           mean +/- 2std = la plage de recherche de ces m/z   (8589 x 5)
    // il ne reste plus qu'a reconstruire les EICs avec ces plages de m/z 
    // --------------------------------------------------------------------------
    

    list_eic4_files=list();
    
    for i=1:n_files;
        
        no_colonnes=mz3.d(:,i,2);
        tri=find(no_colonnes~=0);
        colonnes_gardees=no_colonnes(tri); 
        
        x_tr_mz3=list_eic_files(i); 
        
        x_tr_mz4=x_tr_mz3(:,colonnes_gardees);    
        x_tr_mz4.v=mz3_stats.i(tri);
        x_tr_mz4=div(x_tr_mz4); 
        
        list_eic4_files(i)=x_tr_mz4;
    end

    mz4_stats=mz3_stats;

endfunction