📄 cpnni.pas
字号:
{*********************************************}
{ }
{ COMPONENT for MS DOS and MS WINDOWS }
{ }
{ Source code for Turbo Pascal 6.0 and }
{ Turbo Pasacal for Windows 1.0 compilers. }
{ }
{ (c) 1991, Roderic D. M. Page }
{ }
{*********************************************}
(*{$DEFINE NNIDEBUG}*)
{*
Module to calculate nni metric between two binary trees using
Brown and Day's (1984) heuristic algorithm.
Brown, E. K., and W. H. E. Day. 1984. A computationally efficient
approximation to the nearest neighbour interchange metric.
J. Classif. 1:93-124.
function NNI returns Brown and Day's DRA distance for the
arbitrarily rooted trees. Assumes that the trees are binary.
Routine destroys the trees its is supplied with,
so use with care.
3 Oct 1991 Written and debugged (yeah right!).
10 Oct 1991 Bug found. ExtraNodeCount was initially set to
T1.Leaves + 1. However, if T1 was pruned then
extra nodes would start to overwrite part of T1.NL.
Fixed by making TREEOBJ keep track of how many leaves
had been pruned.
*}
{$I cpdir.inc}
unit cpnni;
interface
uses
{$IFDEF NNIDebug}
WinCrt,
{$ENDIF}
cputil,
cpset,
cptree;
function DRA (var T1, T2: TREEOBJ; rt:integer):integer;
implementation
{ Brown and Day's (1984:Table 2) CARD and CONFIG tables }
const
B = 0;
G = 1;
R = 2;
BG = 3;
BR = 4;
GR = 5;
BGR = 6;
GBR = 7;
RBG = 8;
NBGR = 9;
NGBR = 10;
NRBG = 11;
NNNN = 12;
type
CELL = (CARD,CONFIG);
TABLE = array[B..NNNN,B..NNNN,CELL] of byte;
const
COLORTABLE:TABLE =
(
{B} ((0, B),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0)),
{G} ((0, BG),(0, G),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0)),
{R} ((0, BR),(0, GR),(0, R),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0)),
{BG} ((1, BG),(1, BG),(0, RBG),(3, BG),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0)),
{BR} ((1, BR),(0, GBR),(1, BR),(2,NBGR),(3, BR),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0)),
{GR} ((0, BGR),(1, GR),(1, GR),(2,NGBR),(2,NRBG),(3, GR),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0)),
{BGR} ((1, BGR),(2,NRBG),(2,NGBR),(4,NRBG),(4,NGBR),(4, BGR),(6, BGR),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0)),
{GBR} ((2,NRBG),(1, GBR),(2,NBGR),(4,NRBG),(4, GBR),(4,NBGR),(7,NNNN),(6, GBR),(0, 0),(0, 0),(0, 0),(0, 0),(0, 0)),
{RBG} ((2,NGBR),(2,NBGR),(1, RBG),(4, RBG),(4,NGBR),(4,NBGR),(7,NNNN),(7,NNNN),(6, RBG),(0, 0),(0, 0),(0, 0),(0, 0)),
{NBGR}((2,NNNN),(1, GBR),(1, RBG),(4,NNNN),(4,NNNN),(4,NBGR),(7,NNNN),(6, GBR),(6, RBG),(6,NBGR),(0, 0),(0, 0),(0, 0)),
{NGBR}((1, BGR),(2,NNNN),(1, RBG),(4,NNNN),(4,NGBR),(4,NNNN),(6, BGR),(7,NNNN),(6, RBG),(6, RBG),(6,NGBR),(0, 0),(0, 0)),
{NRBG}((1, BGR),(1, GBR),(2,NNNN),(4,NRBG),(4,NNNN),(4,NNNN),(6, BGR),(6, GBR),(7,NNNN),(6, GBR),(6, BGR),(6,NRBG),(0, 0)),
{NNNN}((1, BGR),(1, GBR),(1, RBG),(4,NNNN),(4,NNNN),(4,NNNN),(6, BGR),(6, GBR),(6, RBG),(6,NBGR),(6,NGBR),(6,NRBG),(6,NNNN))
);
type
SUBTREES = array[0..2] of TREEOBJ_PTR;
var
ExtraNodeCount: integer;
NNICount : integer;
LoopCount : integer;
{-----------------------------ColorTree------------------------------------}
{ Color goaltree S into two or three subtrees, returning
midpoint node (There) and a flag (Edge) for whether
the midpoint is an edge (Edge = true) or a node
(Edge = false). }
procedure ColorTree (var S: TREEOBJ; var There:NODE_PTR; var Edge:Boolean);
var
HalfWeight : integer;
{ Return in t the midpoint of tree S,
Edge is true if midpoint is an edge,
false if it is a node. }
procedure Traverse (p: NODE_PTR);
begin
if (p <> NIL) then begin
Traverse (p^.child);
if (there <> NIL) then
Exit;
if not p^.IsLeaf then
if (p^.GetWeight >= HalfWeight) then begin
there := p;
Edge := (p^.GetWeight = HalfWeight);
Exit;
end;
Traverse (p^.sib);
if (there <> NIL) then
Exit;
end;
end;
{ Paint all leaves above node (Above) with Color }
procedure ColorLeaves (Above: NODE_PTR; Color:Byte);
procedure Traverse (p:NODE_PTR);
begin
if (p <> NIL) then begin
if p^.IsLeaf then
p^.SetIndex (Color);
Traverse (p^.child);
if (p <> Above) then
Traverse (p^.sib);
end;
end;
begin
Traverse (Above);
end;
begin
{ Find midpoint of tree S that breaks S into
2 or 3 more or less equal subtrees. }
There := NIL;
HalfWeight := S.TreeLeaves div 2;
if Odd (S.TreeLeaves) then
Inc (HalfWeight);
Traverse (S.Root);
{ Color all leaves Blue }
ColorLeaves (S.Root, B);
if Edge then
{ S split about an edge, color other subtree Green }
ColorLeaves (There, G)
else begin
{ S split about a node, color left subtree Green,
right subtree Red. }
ColorLeaves (There^.child, G);
ColorLeaves (There^.child^.sib, R);
end;
end;
{-----------------------------LabelT---------------------------------------}
{ Label leaves in start tree T using goal tree S's leaf colors,
assumes that trees are from same profile
(would need to modify for cross
profile comparisons. }
procedure LabelT (var T, S:TREEOBJ);
var
i, j: integer;
begin
j := 0;
i := 0;
while (j < T.TreeLeaves) do begin
repeat
inc (i);
until (T.NL[i] <> NIL);
T.NL[i]^.SetIndex (S.NL[i]^.NodeIndex);
Inc (j);
end;
end;
{-----------------------------Count----------------------------------------}
{ Return the NNI required to coalesce leaves on start tree T of
same color, and return the lub of nonblue leaves in Limit. }
function Count (var T: TREEOBJ; var Limit:NODE_PTR):integer;
var
c,i,j: integer;
LSet :CLUSTEROBJ;
{ Assign colors to internal nodes, and
count nni's. }
procedure Traverse (p:NODE_PTR);
var
i,j,tmp:integer;
begin
if (p <> NIL) then begin
Traverse (p^.child);
if not p^.IsLeaf then begin
i := p^.LeftDesc^.NodeIndex;
j := p^.RightDesc^.NodeIndex;
{ Only lower left of table is occupied }
if (i < j) then begin
tmp := i;
i := j;
j := tmp;
end;
c := c + COLORTABLE[i,j, CARD];
p^.SetIndex (COLORTABLE[i,j, CONFIG]);
end;
if (p <> Limit) then
Traverse (p^.sib);
end;
end;
begin
{ Form set of nonblue leaves }
LSet.NullSet;
{ This loop necessary since subtree leaves
are labelled with labels from original tree
and hence NL[i] is the ith leaf in the original
tree, not necessarily T }
j := 0;
i := 0;
while (j < T.TreeLeaves) do begin
repeat
Inc (i);
until (T.NL[i] <> NIL);
if (T.NL[i]^.NodeIndex <> B) then
LSet.AddToSet (T.NL[i]^.NodeLeafNum);
Inc (j);
end;
{ Find lub of nonblue labels }
Limit := T.lub (LSet);
{ Count nni to transform T into R }
c := 0;
Traverse (Limit);
{ Brown and Day (1984:108), last para in section 2.4:
"if cc(b) is GBR, RBG, or (BGR), then c is incremented
by one since a final nni is required to coalesce the blue
labels above and below b."
}
case Limit^.NodeIndex of
GBR, RBG, NBGR : Inc (c);
else begin end;
end;
Count := c;
end;
{-----------------------------ExtraNode------------------------------------}
{ Return pointer to an extra leaf to be
added to subtree. Node has NamePtr field
set to ExtraNodeCount. }
function ExtraNode:NODE_PTR;
var
p: NODE_PTR;
begin
GetMem (p, SizeOf (TREENODE));
if (p <> NIL) then begin
p^.Init;
with p^ do begin
SetLeaf;
SetNamePtr (ExtraNodeCount);
end;
end;
ExtraNode := p;
end;
{-----------------------------AddExtraNodeGR-------------------------------}
{ Add an extra leaf below the root of a
green or red subtree. }
procedure AddExtraNodeGR (var T:TREEOBJ);
var
p, q: NODE_PTR;
begin
new (p, Init);
q := ExtraNode;
p^.child := q;
q^.anc := p;
q^.sib := T.Root;
T.Root^.anc := p;
T.Root := p;
T.Incleaves;
T.IncInternals;
end;
{-----------------------------AddExtraNodeB--------------------------------}
{ Add an extra leaf to a blue subtree below "Sister" }
procedure AddExtraNodeB (var T:TREEOBJ; Sister: NODE_PTR);
var
p, q, r: NODE_PTR;
begin
new (p, Init);
q := ExtraNode;
p^.child := q;
q^.anc := p;
q^.sib := Sister;
p^.sib := Sister^.sib;
Sister^.sib := NIL;
if Sister^.anc = NIL then
T.Root := p
else begin
if (Sister = Sister^.anc^.child) then
Sister^.anc^.child := p
else Sister^.anc^.child^.sib := p;
end;
p^.anc := Sister^.anc;
Sister^.anc := p;
T.NL[q^.NodeNamePtr] := q;
T.IncLeaves;
T.IncInternals;
end;
{-----------------------------DecompGoalTree-------------------------------}
{ Decompose the goaltree into two or three subtrees of the
same color, depending on whether decomposition is at a node
or an edge.
The original tree T IS altered by this procedure,
becoming pruned to the blue subtree. The blue subtree
is a copy of the pruned T, so that T still exists
in memory to be disposed of later.
}
procedure DecompGoalTree (var T:TREEOBJ; var S: SUBTREES; There: NODE_PTR; Edge:Boolean);
var
S0, S1, S2: TREEOBJ_PTR;
Sister, Ancestor, BelowAnc : NODE_PTR;
LeafCount : integer;
begin
LeafCount := 0;
Ancestor := There^.Anc;
BelowAnc := Ancestor^.anc;
if (There = Ancestor^.Child) then
Sister := There^.Sib
else Sister := Ancestor^.Child;
if Edge then begin
{ Split a tree into two subtrees about edge below There, e.g.:
Colors
B B B G G G
+---Sister
|
1 2 3 4 5 6 1 2 3 4 5 6
\ \ \ \ \ / \ \ \ \ \ /
\ \ \ \ * \ \ \ \ *
\ \ \ \ / \ \ \ \ /
\ \ \ @ <--- There \ \ \ @
\ \ \ / ==> \ \ \
\ \ * <--- Ancestor \ \ *
\ \ / \ \ /
\ * <--- BelowAnc \ *
\ / \ /
* <--- T.Root *
Resulting in
B G
1 2 3 x x 4 5 6
\ \ \ . . \ \ /
\ \ * . \ *
\ \ / . \ /
\ * . @
\ / . .
* *
where x is an extra leaf.
}
{ 1: Green subtree }
{ Prune off subtree by ensuring Sister
is Ancestor's child and has no sibling }
if (There = Sister^.sib) then
Sister^.sib := NIL
else Ancestor^.Child := Sister;
new (S1, Init);
S1^.Root := There;
S1^.Root^.anc := NIL;
{ Ensure green subtree root doesn't have a sister taxon }
S1^.Root^.Sib := NIL;
S1^.SetLeaves (S1^.Root^.GetWeight);
S1^.SetInternals(Pred (S1^.TreeLeaves));
LeafCount := LeafCount + S1^.Treeleaves;
{ Add leaf below root }
AddExtraNodeGR (S1^);
{ 2: Blue subtree }
{ Make the extra leaf a sibling of Sister.
Because we've ensured that Sister is
Ancestor's child (see above) this
is straighforward. }
Sister^.sib := ExtraNode;
Sister^.sib^.anc := Sister^.anc;
T.SetLeaves (T.TreeLeaves - LeafCount + 1);
T.SetInternals (Pred (T.TreeLeaves));
T.NL[Sister^.Sib^.NodeNamePtr] := Sister^.sib;
end
else begin
{ Split into 3 subtrees at node There, e.g. :
Colors
B B B G G G R R
+--- Sister
|
5 6 3 4 1 7 2 8 5 6 3 4 1 7 2 8
\ \ \ \ \ / \ / \ \ \ \ \ / \ /
\ \ \ \ * * \ \ \ \ * *
\ \ \ \ / / \ \ \ \ /
\ \ \ * / \ \ \ *
\ \ \ \ / ==> \ \ \
\ \ \ \ / \ \ \
\ \ \ @ <--- There \ \ \ @
\ \ \ / \ \ \ /
\ \ * <--- Ancestor \ \ *
\ \ / \ \ /
\ * \ *
\ / \ /
* <--- T.Root *
Resulting in
B G R
5 6 3 4 1 7 2 8
\ \ /--x \ \ / \ /
\ * \ * *
\ / \ / |
* * |
| x
|
x
}
{ 1: Green subtree }
new (S1,Init);
S1^.Root := There^.Child;
S1^.Root^.anc := NIL;
S1^.SetLeaves (S1^.Root^.GetWeight);
S1^.SetInternals(Pred (S1^.TreeLeaves));
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -