Skip to content

Commit

Permalink
Merge branch 'intensity-npv-backup' into dev
Browse files Browse the repository at this point in the history
  • Loading branch information
danielgeiszler committed May 29, 2024
2 parents 28cd415 + 776ec32 commit 1ac93ac
Show file tree
Hide file tree
Showing 4 changed files with 180 additions and 25 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -159,14 +159,14 @@ private void fitMatchedIonDistribution() throws Exception {
int mutationSite = decoyPep.mutatedResidue;

// Calculate fake PTM mass shift
float massShift = AAMasses.monoisotopic_masses[targetPep.pepSeq.charAt(mutationSite) - 'A'] -
float mutatedMassShift = AAMasses.monoisotopic_masses[targetPep.pepSeq.charAt(mutationSite) - 'A'] -
AAMasses.monoisotopic_masses[decoyPep.pepSeq.charAt(mutationSite) - 'A'];

// Generate site-determining ions
//XXXQQQ

// Get set of reduced ions that match at least one configuration
float[][] reducedIons = getReducedIons(spec, decoyPep, massShift);
float[][] reducedIons = getReducedIons(spec, decoyPep, mutatedMassShift);
float[] reducedMzs = reducedIons[0];
float[] reducedInts = reducedIons[1];

Expand All @@ -177,26 +177,39 @@ private void fitMatchedIonDistribution() throws Exception {
// decoyPep, this.ionTypes, mutationSite, 1);
// Generate target and decoy fragments
ArrayList<Float> sitePepFrags = targetPep.calculatePeptideFragments(this.ionTypes, 1);
ArrayList<Float> decoySitePepFrags = decoyPep.calculatePeptideFragments(this.ionTypes, 1);
ArrayList<Float> decoySitePepFrags = targetPep.calculateComplementaryFragments(this.ionTypes, mutatedMassShift, mutationSite, 1);

// Score over all possible sites
/**
System.out.println("Decoy\t" + decoyPep.pepSeq + "\t" + Arrays.toString(decoyPep.mods));
System.out.println("Target\t" + targetPep.pepSeq + "\t" + Arrays.toString(targetPep.mods));
System.out.println(mutatedMassShift);
for (int k = 0; k < decoyPep.pepSeq.length(); k++) {
decoyPep.addMod(dMass, k);
ArrayList<Float> sitePepFrags = decoyPep.calculatePeptideFragments(this.ionTypes, 2);
float[][] matchedIons = findMatchedIons(sitePepFrags, reducedMzs, reducedInts);
if (k == mutationSite)
continue;
// Add decoy ions
decoyPep.addMod(mutatedMassShift, k);
ArrayList<Float> decoySitePepFrags = decoyPep.calculatePeptideFragmentsBetween(this.ionTypes, k, mutationSite, 1);
System.out.println("Decoy\t" + k + "\t" + mutationSite + "\t" + decoySitePepFrags);
float[][] matchedIons = findMatchedIons(decoySitePepFrags, reducedMzs, reducedInts);
float[] matchedIonIntensities = matchedIons[0];
System.out.println("Decoy\t" + k + "\t" + mutationSite + "\t" + Arrays.toString(matchedIonIntensities));
float[] matchedIonMassErrors = matchedIons[1];
if (k == mutationSite) {
this.matchedIonDist.addIons(matchedIonIntensities, matchedIonMassErrors, false);
} else {
this.matchedIonDist.addIons(matchedIonIntensities, matchedIonMassErrors, true);
}
decoyPep.addMod(-1*dMass, k);
this.matchedIonDist.addIons(matchedIonIntensities, matchedIonMassErrors, true);
decoyPep.addMod(-1*mutatedMassShift, k);
// Add target ions
ArrayList<Float> sitePepFrags = targetPep.calculatePeptideFragmentsBetween(this.ionTypes, k, mutationSite, 1);
System.out.println("Target\t" + k + "\t" + mutationSite + "\t" + sitePepFrags);
matchedIons = findMatchedIons(sitePepFrags, reducedMzs, reducedInts);
matchedIonIntensities = matchedIons[0];
System.out.println("Target\t" + k + "\t" + mutationSite + "\t" + Arrays.toString(matchedIonIntensities));
matchedIonMassErrors = matchedIons[1];
this.matchedIonDist.addIons(matchedIonIntensities, matchedIonMassErrors, false);
}
**/
// Add target peptide values to matched ion histograms

// Add target peptide values to matched ion histograms
float[][] matchedIons = findMatchedIons(sitePepFrags, reducedMzs, reducedInts);
float[] matchedIonIntensities = matchedIons[0];
float[] matchedIonMassErrors = matchedIons[1];
Expand All @@ -207,6 +220,7 @@ private void fitMatchedIonDistribution() throws Exception {
float[] decoyMatchedIonIntensities = decoyMatchedIons[0];
float[] decoyMatchedMassErrors = decoyMatchedIons[1];
this.matchedIonDist.addIons(decoyMatchedIonIntensities, decoyMatchedMassErrors, true);

}
}
}
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -133,12 +133,16 @@ public void addIon(float intensity, float massError, boolean isDecoy) {
if (intensity > 0.0f) { // Only include matched ions, unmatched ions have negative intensity
this.datapoints.addRow(massError, intensity, isDecoy);
} else {
int intensityIndex = (int) (-1 * intensity);
this.pdf[intensityIndex]++;
this.tNegCount++;
}
} else {
if (intensity > 0.0f) { // Only include matched ions, unmatched ions have negative intensity
this.datapoints.addRow(massError, intensity, isDecoy);
} else {
int intensityIndex = (int) (-1 * intensity);
this.pdfDecoy[intensityIndex]++;
this.dNegCount++;
}
}
Expand Down Expand Up @@ -190,7 +194,7 @@ public void calculateIonPosterior() {
// Sort table and find min and max
this.datapoints.sortRowsByProjVal();
// Set up arrays to hold ion counts for targets and decoys and q-values
int arraySize = (int) (100.0 * (Math.abs(this.datapoints.getMaxProjVal() -
int arraySize = (int) (10.0 * (Math.abs(this.datapoints.getMaxProjVal() -
this.datapoints.getMinProjVal()) + 1));
this.projTargetCounts = new int[arraySize];
this.projDecoyCounts = new int[arraySize];
Expand All @@ -208,17 +212,21 @@ public void calculateIonPosterior() {
boolean leftToRight = checkLdaProjectionLeftToRight();

// Fill q-value array with P(ion == true_positive | LDA score)
int nTargets = 0;
int nTargets = 1;
int nDecoys = 1;
double cMin = 10000.0;
double cmax = -10000.0;
if (leftToRight) {
// Fill array
for (int i = 0; i < arraySize; i++) {
nTargets += this.projTargetCounts[i];
nDecoys += this.projDecoyCounts[i];
this.qVals[i] = (double) nDecoys / (double) (nTargets);
nTargets = this.projTargetCounts[i];
nDecoys = this.projDecoyCounts[i];
double q = (double) (nDecoys + 1) / (double) (nDecoys + nTargets + 2);
this.qVals[i] = q;
System.out.println(this.qVals[i]);
}
// Make array monotonic
/**
for (int i = arraySize-1; i >= 0; i--) {
if (this.qVals[i] < cMin) {
cMin = this.qVals[i];
Expand All @@ -230,13 +238,17 @@ public void calculateIonPosterior() {
i++;
}
}
**/
} else {
for (int i = arraySize-1; i >= 0; i--) {
nTargets += this.projTargetCounts[i];
nDecoys += this.projDecoyCounts[i];
this.qVals[i] = (double) nDecoys / (double) (nTargets);
nTargets = this.projTargetCounts[i];
nDecoys = this.projDecoyCounts[i];
double q = (double) (nDecoys + 1) / (double) (nDecoys + nTargets + 2);
this.qVals[i] = q;
System.out.println(this.qVals[i]);
}
// Make array monotonic
/**
for (int i = 0; i < arraySize; i++) {
if (this.qVals[i] < cMin) {
cMin = this.qVals[i];
Expand All @@ -248,6 +260,7 @@ public void calculateIonPosterior() {
i--;
}
}
**/
}
}

Expand Down Expand Up @@ -281,6 +294,12 @@ public void calculateLdaWeights() throws Exception {

public void calculateNegativePredictiveValue() { // todo see if adding 1 here helps
this.negPredictiveValue = (double) (this.tNegCount + 1) / (double) (this.dNegCount);

this.ionPosterior = new double[((int) (100.0 * this.binMult + 1))];
for (int i = 0; i < this.ionPosterior.length; i++) {
this.ionPosterior[i] = (double) ((this.pdf[i]) / (double) (this.pdf[i] + this.pdfDecoy[i]));
System.out.println(this.ionPosterior[i]);
}
//this.negPredictiveValue = Math.max(this.negPredictiveValue, 0.01);
}

Expand Down Expand Up @@ -388,7 +407,7 @@ public double calcIonProbability(float intensity) {
}

private int translateLdaValToIndex(double projVal) {
int valIndex = (int) ((projVal - this.datapoints.getMinProjVal()) * 100.0);
int valIndex = (int) ((projVal - this.datapoints.getMinProjVal()) * 10.0);
if (valIndex < 0)
return 0;
else if (valIndex >= this.qVals.length)
Expand Down Expand Up @@ -423,9 +442,10 @@ public double calcIonProbability(float intensity, float massError) {
public double calculateIonProbabilityLda(float intensity, float massError) {
double projVal;
double prob;
if (intensity < 0)
if (intensity < 0) {
//prob = this.ionPosterior[(int) (intensity * -1 * this.binMult)];
prob = this.negPredictiveValue;
else {
} else {
projVal = this.ldaProcessor.projectData(massError, intensity).getEntry(0,0);
int projValIndex = translateLdaValToIndex(projVal);
prob = 1.0 - this.qVals[projValIndex];
Expand Down
66 changes: 65 additions & 1 deletion src/edu/umich/andykong/ptmshepherd/utils/Peptide.java
Original file line number Diff line number Diff line change
Expand Up @@ -131,7 +131,71 @@ else if (curIonType == 'x' || curIonType == 'y' || curIonType == 'z')
float cmass = (cTermMass + ccharge * AAMasses.protMass) / ccharge;
for (int i = 0; i < cLen - 1; i++) {
cmass += (aaMasses[this.pepSeq.charAt(cLen - 1 - i) - 'A'] + mods[cLen - 1 - i]) / ccharge;
if (((rightI - 1) >= i) && ( i > (leftI - 1))) {
if ((leftI < (cLen - 1 - i)) && ((cLen - 1 - i) <= rightI)) {
knownFrags.add(cmass);
}
}
}
}

return knownFrags;
}

/**
* Generate peptide ions that are complementary to a modification site (e.g., unshifted where they should be shifted).
* This is designed to be used on a mono-mutated decoy that does not have the mutated mass differene/pseudo delta-mass
* in its modification list.
* @param ionTypes
* @param dmass
* @param site
* @param maxCharge
* @return
*/
public ArrayList<Float> calculateComplementaryFragments(String ionTypes, float dmass, int site, int maxCharge) {
ArrayList<Float> knownFrags = new ArrayList<>(this.pepSeq.length() * ionTypes.length());

// Set up ion types to be identified
ArrayList<Character> nIonTypes = new ArrayList<>();
ArrayList<Character> cIonTypes = new ArrayList<>();
for (int i = 0; i < ionTypes.length(); i++) {
char curIonType = ionTypes.charAt(i);
if (curIonType == 'a' || curIonType == 'b' || curIonType == 'c')
nIonTypes.add(curIonType);
else if (curIonType == 'x' || curIonType == 'y' || curIonType == 'z')
cIonTypes.add(curIonType);
}

float [] aaMasses = AAMasses.monoisotopic_masses;
float [] fragTypeShifts = AAMasses.ionTypeShifts;
int cLen = this.pepSeq.length();

float nTermMass;
for (Character iType : nIonTypes) {
nTermMass = fragTypeShifts[iType - 'a'];
for (int ccharge = 1; ccharge <= maxCharge; ccharge++) { //loop through charge states
float cmass = AAMasses.protMass + nTermMass; //todo this doesn't account for charge?
for (int i = 0; i < cLen - 1; i++) { //loop through positions on the peptide
cmass += (aaMasses[this.pepSeq.charAt(i) - 'A'] + mods[i]) / ccharge;
if (i != 0) { // todo skip a1/b1 ion, check if this needs to be skipped for c too
if (i < site) {
knownFrags.add(cmass + (dmass / ccharge));
} else {
knownFrags.add(cmass);
}
}
}
}
}
float cTermMass;
for (Character iType : cIonTypes) {
cTermMass = fragTypeShifts[iType - 'x' + 3];
for (int ccharge = 1; ccharge <= maxCharge; ccharge++) {
float cmass = (cTermMass + ccharge * AAMasses.protMass) / ccharge;
for (int i = 0; i < cLen - 1; i++) {
cmass += (aaMasses[this.pepSeq.charAt(cLen - 1 - i) - 'A'] + mods[cLen - 1 - i]) / ccharge;
if ((cLen - 1 - i) > site) {
knownFrags.add(cmass + (dmass / ccharge));
} else {
knownFrags.add(cmass);
}
}
Expand Down
57 changes: 57 additions & 0 deletions test/utils/PeptideTest.java
Original file line number Diff line number Diff line change
Expand Up @@ -76,5 +76,62 @@ void calculatePeptideFragmentsBetween() {
for (int i = 0; i < expectedFrags.size(); i++) {
assertEquals(expectedFrags.get(i), sitePepFrags.get(i), 0.0001, "Mismatch at index " + i);
}

expectedFrags = new ArrayList<>(Arrays.asList(227.102633f, 324.155397f));
sitePepFrags = pep.calculatePeptideFragmentsBetween("b", 1, 3, 1);

assertEquals(expectedFrags.size(), sitePepFrags.size());
for (int i = 0; i < expectedFrags.size(); i++) {
assertEquals(expectedFrags.get(i), sitePepFrags.get(i), 0.0001, "Mismatch at index " + i);
}

expectedFrags = new ArrayList<>(Arrays.asList(477.219120f, 574.271884f));
sitePepFrags = pep.calculatePeptideFragmentsBetween("y", 1, 3, 1);

assertEquals(expectedFrags.size(), sitePepFrags.size());
for (int i = 0; i < expectedFrags.size(); i++) {
assertEquals(expectedFrags.get(i), sitePepFrags.get(i), 0.0001, "Mismatch at index " + i);
}
}

@Test
void calculateComplementaryPeptideFragments() {
String seq = "PEPT";
float[] mods = new float[]{0.0f, 0.0f, 0.0f, 0.0f, 0.0f, 0.0f, 0.0f};
Peptide pep = new Peptide(seq, mods);

// Simulate 1.0f dmass on position P3
ArrayList<Float> expectedFrags = new ArrayList<>(Arrays.asList(228.102633f, 324.155397f));
ArrayList<Float> sitePepFrags = pep.calculateComplementaryFragments("b", 1.0f, 2, 1);

assertEquals(expectedFrags.size(), sitePepFrags.size());
for (int i = 0; i < expectedFrags.size(); i++) {
assertEquals(expectedFrags.get(i), sitePepFrags.get(i), 0.0001, "Mismatch at index " + i);
}

// Simulate 1.0f dmass on position P3
expectedFrags = new ArrayList<>(Arrays.asList(121.065520f, 217.118284f, 346.160877f));
sitePepFrags = pep.calculateComplementaryFragments("y", 1.0f, 2, 1);

assertEquals(expectedFrags.size(), sitePepFrags.size());
for (int i = 0; i < expectedFrags.size(); i++) {
assertEquals(expectedFrags.get(i), sitePepFrags.get(i), 0.0001, "Mismatch at index " + i);
}

}

@Test
void AddModTest() {
String seq = "PEP";
float[] mods = new float[]{1.0f, 0.0f, 0.0f};
float[] expectedMods = new float[]{0.0f, 0.0f, 0.0f};
Peptide pep = new Peptide(seq, mods);

pep.addMod(-1*1.0f, 0);
assertArrayEquals(pep.mods, expectedMods);

pep.addMod(1.0f, 0);
expectedMods = new float[]{1.0f, 0.0f, 0.0f};
assertArrayEquals(pep.mods, expectedMods);
}
}

0 comments on commit 1ac93ac

Please sign in to comment.