[Bio] / FigKernelPackages / gjonewicklib.pm Repository:
ViewVC logotype

Diff of /FigKernelPackages/gjonewicklib.pm

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.4, Fri Dec 15 00:09:52 2006 UTC revision 1.5, Tue Dec 19 17:32:48 2006 UTC
# Line 156  Line 156 
156  #  #
157  #  $n_changed = newick_set_undefined_branches( $node, $x )  #  $n_changed = newick_set_undefined_branches( $node, $x )
158  #  $n_changed = newick_set_all_branches( $node, $x )  #  $n_changed = newick_set_all_branches( $node, $x )
159    #  $n_changed = newick_fix_negative_branches( $tree )
160  #  #
161  #  Modify rooting and/or order:  #  Modify rooting and/or order:
162  #  #
# Line 170  Line 171 
171  #  $newtree = reroot_newick_to_node( $tree, @node )  #  $newtree = reroot_newick_to_node( $tree, @node )
172  #  $newtree = reroot_newick_to_node_ref( $tree, $noderef )  #  $newtree = reroot_newick_to_node_ref( $tree, $noderef )
173  #  $newtree = reroot_newick_to_node_ref( $tree, $noderef )  #  $newtree = reroot_newick_to_node_ref( $tree, $noderef )
174    #  $newtree = reroot_newick_to_approx_midpoint( $tree )    # unweighted
175    #  $newtree = reroot_newick_to_approx_midpoint_w( $tree )  # weight by tips
176  #  $newtree = uproot_tip_rooted_newick( $tree )  #  $newtree = uproot_tip_rooted_newick( $tree )
177  #  $newtree = uproot_newick( $tree )  #  $newtree = uproot_newick( $tree )
178  #  #
# Line 177  Line 180 
180  #  $newtree = newick_subtree( $tree,  @tips )  #  $newtree = newick_subtree( $tree,  @tips )
181  #  $newtree = newick_subtree( $tree, \@tips )  #  $newtree = newick_subtree( $tree, \@tips )
182  #  #
183    #  $newtree = collapse_zero_length_branches( $tree )
184    #
185  #===============================================================================  #===============================================================================
186  #  Tree reading and writing:  #  Tree reading and writing:
187  #===============================================================================  #===============================================================================
# Line 235  Line 240 
240    
241          newick_set_undefined_branches          newick_set_undefined_branches
242          newick_set_all_branches          newick_set_all_branches
243            newick_fix_negative_branches
244    
245          normalize_newick_tree          normalize_newick_tree
246          reverse_newick_tree          reverse_newick_tree
# Line 249  Line 255 
255          reroot_newick_to_node          reroot_newick_to_node
256          reroot_newick_to_node_ref          reroot_newick_to_node_ref
257          reroot_newick_to_approx_midpoint          reroot_newick_to_approx_midpoint
258            reroot_newick_to_approx_midpoint_w
259          uproot_tip_rooted_newick          uproot_tip_rooted_newick
260          uproot_newick          uproot_newick
261    
262          prune_from_newick          prune_from_newick
263          newick_subtree          newick_subtree
264            collapse_zero_length_branches
265    
266          writeNewickTree          writeNewickTree
267          fwriteNewickTree          fwriteNewickTree
# Line 1156  Line 1164 
1164    
1165    
1166  #-------------------------------------------------------------------------------  #-------------------------------------------------------------------------------
1167    #  Set negative branches to zero.  The original tree is modfied.
1168    #
1169    #  $n_changed = newick_fix_negative_branches( $tree )
1170    #-------------------------------------------------------------------------------
1171    sub newick_fix_negative_branches {
1172        my ( $tree ) = @_;
1173        array_ref( $tree ) or return undef;
1174        my $n_changed = 0;
1175        my $x = newick_x( $tree );
1176        if ( defined( $x ) and $x < 0 )
1177        {
1178            set_newick_x( $tree, 0 );
1179            $n_changed++;
1180        }
1181    
1182        foreach ( newick_desc_list( $tree ) )
1183        {
1184            $n_changed += newick_fix_negative_branches( $_ );
1185        }
1186    
1187        $n_changed;
1188    }
1189    
1190    
1191    #-------------------------------------------------------------------------------
1192  #  Normalize tree order (in place).  #  Normalize tree order (in place).
1193  #  #
1194  #  ( $tree, $label1 ) = normalize_newick_tree( $tree )  #  ( $tree, $label1 ) = normalize_newick_tree( $tree )
# Line 1500  Line 1533 
1533  #-------------------------------------------------------------------------------  #-------------------------------------------------------------------------------
1534  sub reroot_newick_to_approx_midpoint {  sub reroot_newick_to_approx_midpoint {
1535      my ( $tree ) = @_;      my ( $tree ) = @_;
1536    
1537      #  Compile average tip to node distances assending      #  Compile average tip to node distances assending
1538    
1539      my $dists1 = average_to_tips_1( $tree );      my $dists1 = average_to_tips_1( $tree );
# Line 1563  Line 1597 
1597    
1598      #  The root must be some where below this node:      #  The root must be some where below this node:
1599    
1600      my $n_1      =   @$desc_list + ( $anc_node ? 1 : 0 ) - 1;      my $n_1      =   @$desc_list - ( $anc_node ? 0 : 1 );
1601      my $ttl_dist = ( @$desc_list * $x_below ) + ( defined( $x_above ) ? ( $x_above + $x ) : 0 );      my $ttl_dist = ( @$desc_list * $x_below ) + ( defined( $x_above ) ? ( $x_above + $x ) : 0 );
1602    
1603      foreach ( @$desc_list )      foreach ( @$desc_list )
# Line 1582  Line 1616 
1616    
1617    
1618  #-------------------------------------------------------------------------------  #-------------------------------------------------------------------------------
1619    #  Move root of tree to an approximate midpoint.  Weight by tips.
1620    #
1621    #  $newtree = reroot_newick_to_approx_midpoint_w( $tree )
1622    #-------------------------------------------------------------------------------
1623    sub reroot_newick_to_approx_midpoint_w {
1624        my ( $tree ) = @_;
1625    
1626        #  Compile average tip to node distances assending
1627    
1628        my $dists1 = average_to_tips_1_w( $tree );
1629    
1630        #  Compile average tip to node distances descending, returning midpoint node
1631    
1632        my $node = average_to_tips_2_w( $dists1, undef, undef, undef );
1633    
1634        #  Reroot
1635    
1636        $node ? reroot_newick_to_node_ref( $tree, $node ) : $tree
1637    }
1638    
1639    
1640    sub average_to_tips_1_w {
1641        my ( $node ) = @_;
1642    
1643        my @desc_dists = map { average_to_tips_1_w( $_ ) } newick_desc_list( $node );
1644        my $x_below = 0;
1645        my $n_below = 1;
1646        if ( @desc_dists )
1647        {
1648            $n_below = 0;
1649            my $n;
1650            foreach ( @desc_dists )
1651            {
1652                $n_below += $n = $_->[1];
1653                $x_below += $n * $_->[0];
1654            }
1655            $x_below /= $n_below;
1656        }
1657        my $x = newick_x( $node ) || 0;
1658        my $x_net = $x_below + $x;
1659    
1660        [ $x_net, $n_below, $x, $x_below, [ @desc_dists ], $node ]
1661    }
1662    
1663    
1664    sub average_to_tips_2_w {
1665        my ( $dists1, $x_above, $n_above, $anc_node ) = @_;
1666        my ( undef, $n_below, $x, $x_below, $desc_list, $node ) = @$dists1;
1667    
1668        #  Are we done?  Root is in this node's branch, or "above"?
1669    
1670        # defined( $x_above ) and print STDERR "x_above = $x_above\n";
1671        # print STDERR "x       = $x\n";
1672        # print STDERR "x_below = $x_below\n";
1673        # print STDERR "n_desc  = ", scalar @$desc_list, "\n\n";
1674    
1675        if ( defined( $x_above ) && ( ( $x_above + $x ) >= $x_below ) )
1676        {
1677            #  At this point the root can only be in this node's branch,
1678            #  or "above" it in the current rooting of the tree (which
1679            #  would mean that the midpoint is actually down a different
1680            #  path from the root of the current tree).
1681            #
1682            #  Is the root in the current branch?
1683    
1684            if ( ( $x_below + $x ) >= $x_above )
1685            {
1686                return ( $x_above >= $x_below ) ? $anc_node : $node;
1687            }
1688            else
1689            {
1690                return undef;
1691            }
1692        }
1693    
1694        #  The root must be some where below this node:
1695    
1696        $n_above ||= 0;
1697        my $n = $n_above + $n_below;
1698        my $ttl_w_dist = ( $n_below * $x_below )
1699                       + ( defined( $x_above ) ? $n_above * ( $x_above + $x ) : 0 );
1700    
1701        foreach ( @$desc_list )
1702        {
1703            my $n_2      = $_->[1];    # n in subtree
1704            my $n_above2 = $n - $n_2;  # tip rooted has 1 above
1705    
1706            #  If input tree is tip_rooted, $n_above2 can be 0, so:
1707    
1708            my $x_above2 = $n_above2 ? ( ( $ttl_w_dist - $n_2 * $_->[0] ) / $n_above2 )
1709                                     : 0;
1710            my $root = average_to_tips_2_w( $_, $x_above2, $n_above2 || 1, $node );
1711            if ( $root ) { return $root }
1712        }
1713    
1714        #  Was not anywhere below this node (oh-oh):
1715    
1716        return undef;
1717    }
1718    
1719    
1720    #-------------------------------------------------------------------------------
1721  #  Move root of tree from tip to adjacent node.  #  Move root of tree from tip to adjacent node.
1722  #  #
1723  #  $newtree = uproot_tip_rooted_newick( $tree )  #  $newtree = uproot_tip_rooted_newick( $tree )
# Line 1684  Line 1820 
1820    
1821    
1822  #-------------------------------------------------------------------------------  #-------------------------------------------------------------------------------
1823    #  Collapse zero-length branches to make multifurcation.  The original tree
1824    #  is modified.
1825    #
1826    #  $tree = collapse_zero_length_branches( $tree )
1827    #  $tree = collapse_zero_length_branches( $tree, $not_root )
1828    #-------------------------------------------------------------------------------
1829    sub collapse_zero_length_branches {
1830        my ( $tree, $not_root ) = @_;
1831        array_ref( $tree ) || return undef;
1832    
1833        my @desc = newick_desc_list( $tree );
1834        @desc or return ( $tree );              # Cannot collapse terminal branch
1835    
1836        #  Analyze descendants:
1837    
1838        $not_root ||= 0;
1839        my @new_desc = ();
1840        my $changed = 0;
1841        foreach ( @desc )
1842        {
1843            my ( undef, @to_add ) = collapse_zero_length_branches( $_, $not_root+1 );
1844            if ( @to_add )
1845            {
1846                push @new_desc, @to_add;
1847                $changed = 1;
1848            }
1849            else
1850            {
1851                push @new_desc, $_;
1852            }
1853        }
1854        set_newick_desc_ref( $tree, [ @new_desc ] ) if $changed;
1855    
1856        #  Collapse if not root, not tip and zero (or negative) branch:
1857    
1858        my $collapse = $not_root && @new_desc && ( newick_x( $tree ) <= 0 ) ? 1 : 0;
1859        ( $tree, ( $collapse ? @new_desc : () ) );
1860    }
1861    
1862    
1863    #-------------------------------------------------------------------------------
1864  #  Prune one or more tips from a tree:  #  Prune one or more tips from a tree:
1865  #     Caveat:  if one tip is listed, the original tree is modified.  #     Caveat:  if one tip is listed, the original tree is modified.
1866  #              if more than one tip is listed, a copy of the tree is returen  #              if more than one tip is listed, a copy of the tree is returned
1867  #                   (even if it is just listing the same tip twice!).  #                   (even if it is just listing the same tip twice!).
1868  #  #
1869  #  $newtree = prune_from_newick( $tree,  $tip  )  #  $newtree = prune_from_newick( $tree,  $tip  )
# Line 2205  Line 2382 
2382      if ( ! defined( $file ) ) {      if ( ! defined( $file ) ) {
2383          $fh = \*STDOUT;          $fh = \*STDOUT;
2384      }      }
2385      elsif ( ref($file) eq "GLOB") {      elsif ( $file =~ /^\*/ ) {
2386          $fh = $file;          $fh = $file;
2387      }      }
2388      else {      else {
# Line 2237  Line 2414 
2414    
2415      $min_dx = int( $min_dx );      $min_dx = int( $min_dx );
2416      $dy     = int( $dy );      $dy     = int( $dy );
2417      # RAE: Had an error:      my $x_scale = $width / ( newick_max_X( $node ) || 1 );  # Div by zero caught by RAE
     # Illegal division by zero at /disks/space0/fig/FIGdisk.dev/dist/releases/ross/linux-gentoo/lib/FigKernelPackages/gjonewicklib.pm line 2240, <TREE> line 4., referer: http://listeria.uchicago.edu/dev/FIG/protein.cgi?prot=fig%7C220341.1.peg.3792&sims=1&fid=fig%7C220341.1.peg.3792&user=&translate=0&maxN=500&max_expand=5&maxP=1e-05&select=figx&sort_by=bits&group_by_genome=1&resubmit=resubmit  
     # therefore set newick_max_X to 1 if it is not defined  
   
     my $x_scale = $width / (newick_max_X( $node ) or 1);  
     #my $x_scale = $width / newick_max_X( $node );  
2418    
2419      my $hash = {};      my $hash = {};
2420      layout_printer_plot( $node, $hash, 0, -0.5 * $dy, $x_scale, $min_dx, $dy );      layout_printer_plot( $node, $hash, 0, -0.5 * $dy, $x_scale, $min_dx, $dy );

Legend:
Removed from v.1.4  
changed lines
  Added in v.1.5

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3