Returns standard amino acid BLAST scoring

@param matrix     BLAST substitution matrix
@param gapPenalty penalty for gap, must be < 0
@return standard amino acid BLAST scoring
515234
public static LinearGapAlignmentScoring<AminoAcidSequence> getAminoAcidBLASTScoring(BLASTMatrix matrix, int gapPenalty) {
        return new LinearGapAlignmentScoring<>(AminoAcidSequence.ALPHABET,
                matrix.getMatrix(),
                gapPenalty);
    }
274869
def M(self, t, tips=None, gaps=None):
        """"""
        assert isinstance(t, float) and t > 0, "Invalid t: {0}".format(t)
        with scipy.errstate(under='ignore'): # don't worry if some values 0
            if ('expD', t) not in self._cached:
                self._cached[('expD', t)] = scipy.exp(self.D * self.mu * t)
            expD = self._cached[('expD', t)]
            # swap axes to broadcast multiply D as diagonal matrix
            temp = scipy.ascontiguousarray((self.A.swapaxes(0, 1)
                    * expD).swapaxes(1, 0), dtype=float)
            M = broadcastMatrixMultiply(temp, self.Ainv)
            assert M.min() > -1e-3, "Overly negative M: {0}".format(M.min())
            M[M < 0] = 0.0
            if tips is None:
                return scipy.tile(M, (self.nsites, 1, 1))
            else:
                newM = scipy.zeros((len(tips), N_CODON))
                for i in range(len(tips)):
                    newM[i] =(M[0][:,tips[i]])
                if gaps is not None:
                    newM[gaps] = scipy.ones(N_CODON, dtype='float')
                return newM
716800
private static double getRadiusForAmino(AminoAcid amino, Atom atom) {

		if (atom.getElement().equals(Element.H)) return Element.H.getVDWRadius();
		// some unusual entries (e.g. 1tes) contain Deuterium atoms in standard aminoacids
		if (atom.getElement().equals(Element.D)) return Element.D.getVDWRadius();

		String atomCode = atom.getName();
		char aa = amino.getAminoType();

		// here we use the values that Chothia gives in his paper (as NACCESS does)
		if (atom.getElement()==Element.O) {
			return OXIGEN_VDW;
		}
		else if (atom.getElement()==Element.S) {
			return SULFUR_VDW;
		}
		else if (atom.getElement()==Element.N) {
			if (atomCode.equals("NZ")) return TETRAHEDRAL_NITROGEN_VDW; // tetrahedral Nitrogen
			return TRIGONAL_NITROGEN_VDW;								// trigonal Nitrogen
		}
		else if (atom.getElement()==Element.C) { // it must be a carbon
			if (atomCode.equals("C") ||
					atomCode.equals("CE1") || atomCode.equals("CE2") || atomCode.equals("CE3") ||
					atomCode.equals("CH2") ||
					atomCode.equals("CZ") || atomCode.equals("CZ2") || atomCode.equals("CZ3")) {
				return TRIGONAL_CARBON_VDW; 							// trigonal Carbon
			}
			else if (atomCode.equals("CA") || atomCode.equals("CB") ||
					atomCode.equals("CE") ||
					atomCode.equals("CG1") || atomCode.equals("CG2")) {
				return TETRAHEDRAL_CARBON_VDW;							// tetrahedral Carbon
			}
			// the rest of the cases (CD, CD1, CD2, CG) depend on amino acid
			else {
				switch (aa) {
				case 'F':
				case 'W':
				case 'Y':
				case 'H':
				case 'D':
				case 'N':
					return TRIGONAL_CARBON_VDW;

				case 'P':
				case 'K':
				case 'R':
				case 'M':
				case 'I':
				case 'L':
					return TETRAHEDRAL_CARBON_VDW;

				case 'Q':
				case 'E':
					if (atomCode.equals("CD")) return TRIGONAL_CARBON_VDW;
					else if (atomCode.equals("CG")) return TETRAHEDRAL_CARBON_VDW;

				default:
					logger.info("Unexpected carbon atom "+atomCode+" for aminoacid "+aa+", assigning its standard vdw radius");
					return Element.C.getVDWRadius();
				}
			}

			// not any of the expected atoms
		} else {
			// non standard aas, (e.g. MSE, LLP) will always have this problem,
			logger.info("Unexpected atom "+atomCode+" for aminoacid "+aa+ " ("+amino.getPDBName()+"), assigning its standard vdw radius");


			return atom.getElement().getVDWRadius();
		}
	}
105129
def absolute_extreme_coefficient_ratio(model):
    
    s_matrix, _, _ = con_helpers.stoichiometry_matrix(
        model.metabolites, model.reactions
    )
    abs_matrix = np.abs(s_matrix)
    return abs_matrix.max(), abs_matrix[abs_matrix > 0].min()
1205759
private static float calculateExtinctionFromPeptide(List<Monomer> monomers) throws IOException, HELM2HandledException {

    if (null == monomers || monomers.isEmpty()) {
      return 0.0f;
    }

    Map<String, Integer> countMap = new HashMap<String, Integer>();
    for (Monomer monomer : monomers) {
      if (aminoAcidMap.containsKey(monomer.getAlternateId())) {
        int count = 1;
        if (countMap.containsKey(monomer.getAlternateId())) {
          count = count + countMap.get(monomer.getAlternateId());
        }
        countMap.put(monomer.getAlternateId(), count);
      }
    }

    float result = 0.0f;
    Set<String> keys = countMap.keySet();
    for (Iterator<String> it = keys.iterator(); it.hasNext();) {
      String key = it.next();
      int count = countMap.get(key);
      float factor = aminoAcidMap.get(key);
      result = result + factor * count;
    }

    return BigDecimal.valueOf(result).floatValue();

  }
428899
public double[][] buildNodeMatrix(double[] nodeSequence) {
        int i, j, k;
        for (i = 0; i < nodeNumber; i++) {
            nodeMatrix[i][0] = nodeSequence[i];
            for (j = 1; j <= layerNumber; j++) {
                nodeMatrix[i][j] = 0.0;
                for (k = 0; k < nodeNumber; k++) {
                    if (apspMatrix[i][k] == j) {
                        nodeMatrix[i][j] += nodeSequence[k];
                    }
                }
            }
        }
        return nodeMatrix;
    }
941686
function calculatePenalty(matrix) {
    var N = matrix.length;
    var penalty = 0;
    // Rule 1
    for (var i = 0; i < N; i++) {
        var pixel = matrix[i][0] & 1;
        var len = 1;
        for (var j = 1; j < N; j++) {
            var p = matrix[i][j] & 1;
            if (p == pixel) {
                len++;
                continue;
            }
            if (len >= 5) {
                penalty += len - 2;
            }
            pixel = p;
            len = 1;
        }
        if (len >= 5) {
            penalty += len - 2;
        }
    }
    for (var j = 0; j < N; j++) {
        var pixel = matrix[0][j] & 1;
        var len = 1;
        for (var i = 1; i < N; i++) {
            var p = matrix[i][j] & 1;
            if (p == pixel) {
                len++;
                continue;
            }
            if (len >= 5) {
                penalty += len - 2;
            }
            pixel = p;
            len = 1;
        }
        if (len >= 5) {
            penalty += len - 2;
        }
    }

    // Rule 2
    for (var i = 0; i < N - 1; i++) {
        for (var j = 0; j < N - 1; j++) {
            var s = matrix[i][j] + matrix[i][j + 1] + matrix[i + 1][j] + matrix[i + 1][j + 1] & 7;
            if (s == 0 || s == 4) {
                penalty += 3;
            }
        }
    }

    // Rule 3
    function I(k) { return matrix[i][j + k] & 1 };
    function J(k) { return matrix[i + k][j] & 1 };
    for (var i = 0; i < N; i++) {
        for (var j = 0; j < N; j++) {
            if (j < N - 6 && I(0) && !I(1) && I(2) && I(3) && I(4) && !I(5) && I(6)) {
                if (j >= 4 && !(I(-4) || I(-3) || I(-2) || I(-1))) {
                    penalty += 40;
                }
                if (j < N - 10 && !(I(7) || I(8) || I(9) || I(10))) {
                    penalty += 40;
                }
            }

            if (i < N - 6 && J(0) && !J(1) && J(2) && J(3) && J(4) && !J(5) && J(6)) {
                if (i >= 4 && !(J(-4) || J(-3) || J(-2) || J(-1))) {
                    penalty += 40;
                }
                if (i < N - 10 && !(J(7) || J(8) || J(9) || J(10))) {
                    penalty += 40;
                }
            }
        }
    }

    // Rule 4
    var numDark = 0;
    for (var i = 0; i < N; i++) {
        for (var j = 0; j < N; j++) {
            if (matrix[i][j] & 1) numDark++;
        }
    }
    penalty += 10 * Math.floor(Math.abs(10 - 20 * numDark/(N * N)));

    return penalty;
}
43520
public static DoubleMatrix cholesky(DoubleMatrix A) {
      DoubleMatrix result = A.dup();
      int info = NativeBlas.dpotrf('U', A.rows, result.data, 0, A.rows);
      if (info < 0) {
          throw new LapackArgumentException("DPOTRF", -info);
      } else if (info > 0) {
          throw new LapackPositivityException("DPOTRF", "Minor " + info + " was negative. Matrix must be positive definite.");
      }
      clearLower(result);
      return result;
  }
564823
public static DoubleMatrix ger(double alpha, DoubleMatrix x,
			DoubleMatrix y, DoubleMatrix a) {
		NativeBlas.dger(a.rows, a.columns, alpha, x.data, 0, 1, y.data, 0, 1, a.data,
				0, a.rows);
		return a;
	}
926198
public int getBlastGapExtension() {
		String gapCosts = getAlignmentOption(GAPCOSTS);
		try {
			String gapExtension = gapCosts.split("\\+")[1];
			return Integer.parseInt(gapExtension);
		} catch (Exception e) {
			return -1;
		}
	}
1209275
function getMatrix(data) {
    var matrix = init(data.version);
    fillFinders(matrix);
    fillAlignAndTiming(matrix);
    fillStub(matrix);

    var penalty = Infinity;
    var bestMask = 0;
    for (var mask = 0; mask < 8; mask++) {
        fillData(matrix, data, mask);
        fillReserved(matrix, data.ec_level, mask);
        var p = calculatePenalty(matrix);
        if (p < penalty) {
            penalty = p;
            bestMask = mask;
        }
    }

    fillData(matrix, data, bestMask);
    fillReserved(matrix, data.ec_level, bestMask);

    return matrix.map(function(row) {
        return row.map(function(cell) {
            return cell & 1;
        });
    });
}
85435
public static <S extends Sequence<S>> int calculateScore(S seq1, Range seq1Range, Mutations<S> mutations,
                                                             LinearGapAlignmentScoring<S> scoring) {
        if (!mutations.isEmpty() && mutations.getPositionByIndex(0) < seq1Range.getFrom() - 1)
            throw new IllegalArgumentException();

        final AlignmentIteratorForward<S> iterator = new AlignmentIteratorForward<>(mutations, seq1Range);

        int score = 0;

        while (iterator.advance()) {
            final int mut = iterator.getCurrentMutation();
            switch (Mutation.getRawTypeCode(mut)) {
                case RAW_MUTATION_TYPE_SUBSTITUTION:
                    score += scoring.getScore(Mutation.getFrom(mut), Mutation.getTo(mut));
                    break;

                case RAW_MUTATION_TYPE_DELETION:
                case RAW_MUTATION_TYPE_INSERTION:
                    score += scoring.getGapPenalty();
                    break;

                default:
                    byte c = seq1.codeAt(iterator.getSeq1Position());
                    score += scoring.getScore(c, c);
                    break;
            }
        }

        return score;
    }
464789
static HashMap<String, BlastResults> readBlast(final int alleleCutoff, final File inputBlastFile) throws IOException{
    	
        BufferedReader reader = null;
        HashMap<String, BlastResults> blast = new HashMap<String, BlastResults>();
        
        try {
            reader = new BufferedReader(new FileReader(inputBlastFile));

            for (HighScoringPair hsp : HspReader.read(reader)) {

            	//id	allele 
            	List<String> BlastnParts = Splitter.on("\t").splitToList(hsp.toString());
            	String sequenceId        = BlastnParts.get(0);
            	String allele            = BlastnParts.get(1);
            
               BlastResults blastResults    = blast.get(sequenceId);
               if(blastResults == null){
            	   blastResults  = new BlastResults(alleleCutoff);
               	   blast.put(sequenceId, blastResults);
               }
               
               blastResults.addTyping(allele);
            }
            
        }finally {
            try {
            	reader.close();
            }
            catch (Exception e) {
            	 e.printStackTrace();
            	 System.exit(1);
            }
        }
    
        
        return blast;
    	
    }
1058525
public static <S extends Sequence<S>> Alignment<S> alignGlobal(AlignmentScoring<S> alignmentScoring,
                                                                   S seq1, S seq2) {
        if (alignmentScoring instanceof AffineGapAlignmentScoring)
            return alignGlobalAffine((AffineGapAlignmentScoring<S>) alignmentScoring, seq1, seq2);
        if (alignmentScoring instanceof LinearGapAlignmentScoring)
            return alignGlobalLinear((LinearGapAlignmentScoring<S>) alignmentScoring, seq1, seq2);
        throw new RuntimeException("Unknown scoring type.");
    }
383119
public int calculateScore(AlignmentScoring<S> scoring) {
        return AlignmentUtils.calculateScore(sequence1, sequence1Range, mutations, scoring);
    }
282120
public String[] getBestLabels(double[] scores, float ratioOfBest) {
		if (ratioOfBest > 1 | ratioOfBest <= 0)
			throw new RuntimeException(
					"ratioOfBest shoudl be > 0 and <= 1 but got " + ratioOfBest);
		TreeMap<Double, String> sortedMap = new TreeMap<Double, String>();
		double max = -Double.MAX_VALUE;
		double min = Double.MAX_VALUE;
		for (int d = 0; d < scores.length; d++) {
			if (scores[d] > max) {
				max = scores[d];
			}
			if (scores[d] < min) {
				min = scores[d];
			}
		}

		double scaleFactor = max - min;
		for (int x = 0; x < scores.length; x++) {
			double newScore = ((scores[x] - min) / scaleFactor);
			sortedMap.put(newScore, this.lexicon.getLabel(x));
		}

		// find cutoffpoint -> we know that the best label will have a value of 1
		double threshold = ratioOfBest;

		List<String> labelsKept = new ArrayList<String>();
		Iterator<Entry<Double, String>> pairs = sortedMap.entrySet().iterator();
		while (pairs.hasNext()) {
			Entry<Double, String> pair = pairs.next();
			if (pair.getKey() >= threshold)
				labelsKept.add(pair.getValue());
		}
		return (String[]) labelsKept.toArray(new String[labelsKept.size()]);
	}
1229679
public static Map.Entry<Object, Object> maxMax(DataTable2D payoffMatrix) {
        if(payoffMatrix.isValid()==false) {
            throw new IllegalArgumentException("The payoff matrix does not have a rectangular format.");
        }
        
        Double maxMaxPayoff = Double.NEGATIVE_INFINITY;
        Object maxMaxPayoffOption = null;
        
        for(Map.Entry<Object, AssociativeArray> entry : payoffMatrix.entrySet()) {
            //Object event = entry.getKey();
            AssociativeArray optionList = entry.getValue();
            for(Map.Entry<Object, Object> entry2 : optionList.entrySet()) {
                Double payoff = TypeInference.toDouble(entry2.getValue());
                //Object option = entry2.getKey();
                
                if(payoff>maxMaxPayoff) {
                    maxMaxPayoff=payoff;
                    maxMaxPayoffOption = entry2.getKey();
                }
            }
        }
        
        return new AbstractMap.SimpleEntry<>(maxMaxPayoffOption, maxMaxPayoff);
    }