X-Git-Url: http://source.jalview.org/gitweb/?a=blobdiff_plain;f=forester%2Fjava%2Fsrc%2Forg%2Fforester%2Fevoinference%2Fdistance%2FNeighborJoining.java;h=bcb6df1e09ebcbafe8b6715dce228d96568e481c;hb=2eb16fa8fc6fa40a6f9829a5d624ada26b026fcc;hp=ef7fe920259a154beb92cbb93261f7794aa7a3d4;hpb=e2785e20355846acd9e81d50cf9a8a56b04ec095;p=jalview.git diff --git a/forester/java/src/org/forester/evoinference/distance/NeighborJoining.java b/forester/java/src/org/forester/evoinference/distance/NeighborJoining.java index ef7fe92..bcb6df1 100644 --- a/forester/java/src/org/forester/evoinference/distance/NeighborJoining.java +++ b/forester/java/src/org/forester/evoinference/distance/NeighborJoining.java @@ -21,7 +21,7 @@ // Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA // // Contact: phylosoft @ gmail . com -// WWW: www.phylosoft.org/forester +// WWW: https://sites.google.com/site/cmzmasek/home/software/forester package org.forester.evoinference.distance; @@ -36,7 +36,9 @@ import org.forester.util.ForesterUtil; public final class NeighborJoining { private BasicSymmetricalDistanceMatrix _d; - private BasicSymmetricalDistanceMatrix _m; + // private BasicSymmetricalDistanceMatrix _m; + private double[][] _d_values; + private double[][] _m_values; private double[] _r; private int _n; private PhylogenyNode[] _external_nodes; @@ -49,24 +51,27 @@ public final class NeighborJoining { } private final void calculateDistancesFromNewNode( final int otu1, final int otu2, final double d ) { + final int otu1_m = _mappings[ otu1 ]; + final int otu2_m = _mappings[ otu2 ]; + int i_m; for( int i = 0; i < _n; ++i ) { if ( ( i == otu1 ) || ( i == otu2 ) ) { continue; } - //final double nd = ; - // setValueInD( nd, otu1, i ); - // _d.setValue( _mappings[ otu1 ], - // _mappings[ i ], - // ( getValueFromD( otu1, i ) + getValueFromD( i, otu2 ) - d ) / 2 ); - _d._values[ _mappings[ otu1 ] ][ _mappings[ i ] ] = ( getValueFromD( otu1, i ) + getValueFromD( i, otu2 ) - d ) / 2; + _d_values[ _mappings[ otu1 ] ][ _mappings[ i ] ] = ( getValueFromD( otu1, i ) + getValueFromD( i, otu2 ) - d ) / 2; + //i_m = _mappings[ i ]; + //_d_values[ otu1_m ][ i_m ] = ( ( _d_values[ otu1_m ][ i_m ] + _d_values[ i_m ][ otu2_m ] ) - 2 ) / 2; } } private final void calculateNetDivergences() { double d; + int i_m; for( int i = 0; i < _n; ++i ) { d = 0; + i_m = _mappings[ i ]; for( int n = 0; n < _n; ++n ) { + //d += _d_values[ i_m ][ _mappings[ n ] ]; d += getValueFromD( i, n ); } _r[ i ] = d; @@ -78,7 +83,6 @@ public final class NeighborJoining { final Phylogeny phylogeny = new Phylogeny(); while ( _n > 2 ) { updateM(); - //final int[] s = findMinimalDistance(); // Calculates the minimal distance. // If more than one minimal distances, always the first found is used // could randomize this, so that any would be returned in a randomized fashion... @@ -87,22 +91,17 @@ public final class NeighborJoining { int otu2 = -1; for( int j = 1; j < _n; ++j ) { for( int i = 0; i < j; ++i ) { - // if ( _m.getValue( i, j ) < minimum ) { - if ( _m._values[ i ][ j ] < minimum ) { - //minimum = _m.getValue( i, j ); - minimum = _m._values[ i ][ j ]; + if ( _m_values[ i ][ j ] < minimum ) { + minimum = _m_values[ i ][ j ]; otu1 = i; otu2 = j; } } } - // - // final int otu1 = s[ 0 ]; - // final int otu2 = s[ 1 ]; // It is a condition that otu1 < otu2. if ( DEBUG ) { if ( otu1 > otu2 ) { - throw new RuntimeException( "NJ code is faulty: otu1 > otu2" ); + throw new RuntimeException( "faulty NJ code: otu1 > otu2" ); } } final PhylogenyNode node = new PhylogenyNode(); @@ -117,7 +116,6 @@ public final class NeighborJoining { printProgress( otu1, otu2 ); } calculateDistancesFromNewNode( otu1, otu2, d ); - //setExternalPhylogenyNode( node, otu1 ); _external_nodes[ _mappings[ otu1 ] ] = node; updateMappings( otu2 ); --_n; @@ -144,39 +142,20 @@ public final class NeighborJoining { return pl; } - // private int[] findMinimalDistance() { - // // if more than one minimal distances, always the first found is - // // returned - // // i could randomize this, so that any would be returned in a randomized - // // fashion... - // double minimum = Double.MAX_VALUE; - // int otu_1 = -1; - // int otu_2 = -1; - // for( int j = 1; j < _n; ++j ) { - // for( int i = 0; i < j; ++i ) { - // if ( _m.getValue( i, j ) < minimum ) { - // minimum = _m.getValue( i, j ); - // otu_1 = i; - // otu_2 = j; - // } - // } - // } - // return new int[] { otu_1, otu_2 }; - // } private final PhylogenyNode getExternalPhylogenyNode( final int i ) { return _external_nodes[ _mappings[ i ] ]; } private final double getValueFromD( final int otu1, final int otu2 ) { - //return _d.getValue( _mappings[ otu1 ], _mappings[ otu2 ] ); - return _d._values[ _mappings[ otu1 ] ][ _mappings[ otu2 ] ]; + return _d_values[ _mappings[ otu1 ] ][ _mappings[ otu2 ] ]; } private final void initExternalNodes() { _external_nodes = new PhylogenyNode[ _n ]; + String id; for( int i = 0; i < _n; ++i ) { _external_nodes[ i ] = new PhylogenyNode(); - final String id = _d.getIdentifier( i ); + id = _d.getIdentifier( i ); if ( id != null ) { _external_nodes[ i ].setName( id ); } @@ -199,35 +178,34 @@ public final class NeighborJoining { private final void reset( final BasicSymmetricalDistanceMatrix distances ) { _n = distances.getSize(); _d = distances; - _m = new BasicSymmetricalDistanceMatrix( _n ); _r = new double[ _n ]; _mappings = new int[ _n ]; + _d_values = _d.getValues(); + _m_values = new double[ _n ][ _n ]; initExternalNodes(); } - // private final void setExternalPhylogenyNode( final PhylogenyNode node, final int i ) { - // _external_nodes[ _mappings[ i ] ] = node; - // } - // private final void setValueInD( final double d, final int otu1, final int otu2 ) { - // _d.setValue( _mappings[ otu1 ], _mappings[ otu2 ], d ); - // } private final void updateM() { calculateNetDivergences(); + double r_j; + int j_m; + final int _n_2 = _n - 2; for( int j = 1; j < _n; ++j ) { + r_j = _r[ j ]; + j_m = _mappings[ j ]; for( int i = 0; i < j; ++i ) { - //_m.setValue( i, j, calculateM( i, j ) ); - //_m.setValue( i, j, getValueFromD( i, j ) - ( _r[ i ] + _r[ j ] ) / ( _n - 2 ) ); - _m._values[ i ][ j ] = getValueFromD( i, j ) - ( _r[ i ] + _r[ j ] ) / ( _n - 2 ); + _m_values[ i ][ j ] = getValueFromD( i, j ) - ( _r[ i ] + r_j ) / ( _n - 2 ); + //_m_values[ i ][ j ] = _d_values[ _mappings[ i ] ][ j_m ] - ( ( _r[ i ] + r_j ) / ( _n_2 ) ); } } } - //private double calculateM( final int i, final int j ) { - // return getValueFromD( i, j ) - ( _r[ i ] + _r[ j ] ) / ( _n - 2 ); - //} + // private final double getValueFromD( final int otu1, final int otu2 ) { + // return _d_values[ _mappings[ otu1 ] ][ _mappings[ otu2 ] ]; + // } // otu2 will, in effect, be "deleted" from the matrix. private final void updateMappings( final int otu2 ) { - for( int i = otu2; i < _mappings.length - 1; ++i ) { + for( int i = otu2; i < ( _mappings.length - 1 ); ++i ) { _mappings[ i ] = _mappings[ i + 1 ]; } }