
(*==============================================================================

     composition.tl

     This file is part of the TREEFINDER software.

                                    COPYRIGHT (C) 1997-2008 BY GANGOLF JOBB
==============================================================================*)

CheckBaseComposition[
 data,
 WindowSize->500,
 WindowStep->100
]:=<<
 Scope[{counts,total,namewidth,composition,specific,names,line},

  data=RemoveFilters[LoadSequences[data]];

  counts=BaseCounts[data];
  If[counts@Mul==0,Throw["CheckBaseComposition: Vanishing nucleotide frequency."]];
  total=counts/counts@Add;
  AppendR[total,total|1+total|2];
  AppendR[total,total|3+total|4];

  namewidth=4+Max[GetSequenceNames[data]#StringLength];
  If[namewidth<8,namewidth=8];

  composition[sequence]=
   <<Scope[{c:=BaseCounts[Flatten[{sequence}]]},c=c/c@Add,Append[Append[c,c|1+c|2],c|3+c|4]]>>;

  specific=Transpose[data#composition];
  names=data#<<1,Get>>,

  line[sequence]=<<
   "  "
   <>PadRight[ToString[sequence|2],namewidth]
   <>":  "
   <>sequence|1
  >>;

  Join[

   {"",
    "BASE COMPOSITION",
    "================",
    "",
    "Data size:",
    "",
    "  "<>Length[data]<>" sequences",
    "  "<>Max[GetSequenceLengths[data]]<>" sites (max.)",
    "",
    "Total base composition:",
    "",
    "  T   : "<>total|1<>" ("<>counts|1<>")",
    "  C   : "<>total|2<>" ("<>counts|2<>")",
    "  A   : "<>total|3<>" ("<>counts|3<>")",
    "  G   : "<>total|4<>" ("<>counts|4<>")",
    "  T+C : "<>total|5<>" ("<>counts|1+counts|2<>")",
    "  A+G : "<>total|6<>" ("<>counts|3+counts|4<>")",
    "",
    "Sequence specific base composition {T,C,A,G,T+C,A+G}:",
    ""
   },

   Transpose[{Transpose[specific],names}]#line,

   {"",
    "Sequences sorted by T content:",
    ""
   },

   Sort[Transpose[{specific|1,names}]]#line,

   {"",
    "Sequences sorted by C content:",
    ""
   },

   Sort[Transpose[{specific|2,names}]]#line,

   {"",
    "Sequences sorted by A content:",
    ""
   },

   Sort[Transpose[{specific|3,names}]]#line,

   {"",
    "Sequences sorted by G content:",
    ""
   },

   Sort[Transpose[{specific|4,names}]]#line,

   {"",
    "Sequences sorted by T+C content:",
    ""
   },

   Sort[Transpose[{specific|5,names}]]#line,

   {"",
    "Sequences sorted by A+G content:",
    ""
   },

   Sort[Transpose[{specific|6,names}]]#line,

   {""
   }

  ];EditLines

  ;
  Scope[{size:=WindowSize,step:=WindowStep,end,indices,parts,profile,addpos},
   end=Min[GetSequenceLengths[data]];
   indices={};
   parts=
    Table[
     CutSequences[data,{i,j}];AppendR[indices,i]
     ,(=i,1:end:step,j,size:end:step=)
    ];
   If[Length[parts]>1,
    profile=Transpose[parts#composition];
    addpos[curve]=<<Transpose[{indices,curve}]>>;
    PlotPoints[
     profile#addpos,
     DiagramHeight->121.,
     XRange->1:Max[indices],
     YRange->0:0.9,
     Joined->True,
     Color->{ 1:0:0 , 0:0:0.8 , 0:0.9:0 , 0.6:0.4:0.2 , 1:0:0.8 , 1:0.9:0 },
     XLabel->"window start",
     YLabel->"base frequency",
     Header->"Nucleotide Frequency Distribution Over Sites",
     Description->"      red=T  blue=C  green=A  brown=G  magenta=T+C  yellow=A+G\n\n"
                <>"      window size: "<>size<>"      window step: "<>step
    ]
   ]
  ]

  ;
  If[Length[data]>1,
   PlotValues[
    specific#Sort,
    DiagramHeight->121.,
    XRange->1:Length[data],
    YRange->0:0.9,
    Joined->True,
    Color->{ 1:0:0 , 0:0:0.8 , 0:0.9:0 , 0.6:0.4:0.2 , 1:0:0.8 , 1:0.9:0 },
    XLabel->"sorted species",
    YLabel->"base frequency",
    Header->"Nucleotide Frequency Distribution Over Species",
    Description->"      red=T  blue=C  green=A  brown=G  magenta=T+C  yellow=A+G\n\n"
   ]
  ]

 ]
>>,
SetMask[CheckBaseComposition,"1."], Protect[CheckBaseComposition],

(*============================================================================*)

CheckNameCompatibility[
 d1,
 d2
]:=<<

 If[
  Catch[d1=LoadTree[d1]]=="",d1=GetTipNames[d1],
  Catch[d1=LoadSequences[d1]]=="",d1=GetSequenceNames[d1],
  Throw["CheckNameCompatibility: Bad data."]
 ],

 If[
  Catch[d2=LoadTree[d2]]=="",d2=GetTipNames[d2],
  Catch[d2=LoadSequences[d2]]=="",d2=GetSequenceNames[d2],
  Throw["CheckNameCompatibility: Bad data."]
 ],

 Join[
  {"",
   "NAME COMPATIBILITY",
   "==================",
   "",
   "Names that appear in the first\nfile but not in the second:",
   ""
  },
  d1^|{_a=>If[!IsMember[d2,_a],_a]},'ToString,Map1,".",Append,
  {"",
   "Names that appear in the second\nfile but not in the first:",
   ""
  },
  d2^|{_a=>If[!IsMember[d1,_a],_a]},'ToString,Map1,".",Append,
  {""
  }
 ],EditLines

>>,
SetMask[CheckNameCompatibility,"11"], Protect[CheckNameCompatibility],

(*============================================================================*)

()

