image thumbnail

Ancestral polytree

by

 

01 Feb 2013 (Updated )

The code here provides the algorithm of learning ancestral polytree.

LAP_Example.m
%function LAP_Example()
clc
s = matlabpool('size');
if s == 0, matlabpool open; end

cd ./
mex ComputeMutualInformationMatrix.c
mex UndirectedMaximumSpanningTree.c



inputfile =  '/data/HIVB_Capsid.fasta';
algn = ExtractSequenceOut( inputfile );
TransferInt = TransferAA2Int( algn,1 );
MIValue  = ComputeMutualInformationMatrix( TransferInt );

Cutoff = 0.015;
for p = 1:size(MIValue,1)
    MIValue(p,p) = 0;
    for q = 1:size(MIValue,2)
        if MIValue(p,q) < Cutoff
           MIValue(p,q) = 0;
        end
    end
end

PolytreeSkeleton = UndirectedMaximumSpanningTree( MIValue );
% h4 = view(biograph( PolytreeSkeleton )) 
Dim = size(PolytreeSkeleton);
if 1
SelectVariable = zeros( 1,size(PolytreeSkeleton) );
for p = 1:length(PolytreeSkeleton)
    PolytreeSkeleton(p,p) = 0;
    if sum( PolytreeSkeleton(p,:) ) > 0
       SelectVariable(p)=1;
    end
end
SelectVariable = find(SelectVariable==1);
end

Cluster ={1:15,[41 120 132 136 91 86 92  87 96 98 116]};
for p = 2:length(Cluster)
 %   SelectVariable = Cluster{p};
   LGObj = ConstructLGObj( TransferInt(:,SelectVariable) ); 
  % keep the significance level at 0.05;
  CIMatric = ConditionalIndependentMatric( LGObj,0.05 );
  [ Max_UPRMatric,Max_IAPS,Max_CINumber,Max_AmbiguousArc,Success1 ]  = LAP_PolytreeRecover( CIMatric,PolytreeSkeleton(SelectVariable,SelectVariable)>0);
  OrigMatrix = zeros( Dim );
  OrigMatrix(SelectVariable,SelectVariable) = Max_UPRMatric;
  h4 = view(biograph( OrigMatrix )) 
d
end



Contact us