Example program converting ms_ms2quantitation to ms_customquantitation.
import msparser
import sys
if len(sys.argv) != 5:
print("Usage: %s <quantitation schema 1 filepath> <quantitation schema 2 filepath> <unimod schema filepath> <results file>" % sys.argv[0])
print("The location of schema files should be e.g.")
print(" ../html/xmlns/schema/quantitation_1/quantitation_1.xsd")
print(" ../html/xmlns/schema/quantitation_2/quantitation_2.xsd")
print(" ../html/xmlns/schema/unimod_2/unimod_2.xsd")
print("if running from the Mascot cgi directory.")
sys.exit(1)
quant_schema1_filepath = sys.argv[1]
quant_schema2_filepath = sys.argv[2]
unimod_schema_filepath = sys.argv[3]
resfile = msparser.ms_mascotresfile(sys.argv[4], 1)
if not resfile.isValid():
print(resfile.getLastErrorString())
sys.exit(1)
def ms_range(start, stop, step=1):
i = start
while i <= stop:
yield i
i += step
def dump_warnings(err):
for i in ms_range(1, err.getNumberOfErrors()):
print('Warning: %s' %(err.getErrorString(i)))
dump_warnings(resfile.getErrorHandler())
resfile.clearAllErrors()
resfile.setXMLschemaFilePath(msparser.ms_mascotresfile.XML_SCHEMA_QUANTITATION, "http://www.matrixscience.com/xmlns/schema/quantitation_2 " + quant_schema2_filepath + " http://www.matrixscience.com/xmlns/schema/quantitation_1 " + quant_schema1_filepath)
resfile.setXMLschemaFilePath(msparser.ms_mascotresfile.XML_SCHEMA_UNIMOD, "http://www.unimod.org/xmlns/schema/unimod_2 " + unimod_schema_filepath)
qmethod = msparser.ms_quant_method()
if not resfile.getQuantitationMethod(qmethod):
print(resfile.getLastErrorString())
sys.exit(1)
datfile = msparser.ms_datfile()
undef, flags, minprob, maxhits, iisb, minpeplen, use_pepsum, flags2 = resfile.get_ms_mascotresults_params(datfile.getMascotOptions())
if not use_pepsum:
print("Results File cannot be opened as a peptide summary")
sys.exit(1)
pepsum = msparser.ms_peptidesummary(resfile, flags, minprob, maxhits, '', iisb, minpeplen, '', flags2)
if not resfile.isValid():
print(resfile.getLastErrorString)
sys.exit(1)
dump_warnings(resfile.getErrorHandler())
resfile.clearAllErrors()
ms2quant = msparser.ms_ms2quantitation(pepsum, qmethod)
if not ms2quant.isValid():
print(ms2quant.getLastErrorString())
sys.exit(1)
if ms2quant.hasIntensityNormalisation():
ms2quant.normaliseIntensities()
elif ms2quant.hasPeptideRatioNormalisation():
ms2quant.normalisePeptideRatios()
dump_warnings(ms2quant.getErrorHandler())
ms2quant.clearAllErrors()
customquant = msparser.ms_customquantitation(ms2quant)
if not customquant.isValid():
print(customquant.getLastErrorString())
sys.exit(1)
err = customquant.getErrorHandler()
for i in range(1, err.getNumberOfErrors()):
print("Warning: %s" %(err.getErrorString(i)))
customquant.clearAllErrors()
qmethod = customquant.getQuantitationMethod()
comps = []
for i in ms_range(0, qmethod.getNumberOfComponents()-1):
comp = qmethod.getComponentByNumber(i)
comps.append(comp.getName())
print("Components: %s" %(comps))
rationames = []
for i in ms_range(0, qmethod.getNumberOfReportRatios()-1):
ratioByNumber = qmethod.getReportRatioByNumber(i)
rationames.append(ratioByNumber.getName())
print("Report Ratios: %s" %(rationames))
print("Protein ratio type = %s" %(qmethod.getProteinRatioType()))
print("Min. number of peptides = %d" %(qmethod.getMinNumPeptides()))
if qmethod.haveQuality():
q = qmethod.getQuality()
print("Quality: min. precursor charge = %s" %(q.getMinPrecursorCharge()))
print("Quality: pep. threshold type = %s" %(q.getPepThresholdType()))
print("Quality: pep. threshold value = %s" %(q.getPepThresholdValue()))
else:
print("Quality: no restrictions")
if qmethod.haveNormalisation():
q = qmethod.getNormalisation()
print("Normalisation = %s" %(q.getMethod()))
if q.haveProteins():
p = q.getProteins()
print("getProteins")
for i in ms_range(0, p.getNumberOfProteins()-1):
proteins = p.getProteins(i)
protein = proteins.getAccession()
print(" - Accessions = %s" %(protein))
if q.havePeptides():
p = q.getPeptides()
peptides = p.getPeptides()
for i in ms_range(0, p.getNumberOfPeptides()-1):
peptides = p.getPeptides(i)
peptide = peptides.getSequence()
print(" - Sequences = %s" %(peptide))
component_bases = []
for i in ms_range(0, qmethod.getNumberOfComponents()-1):
gcbn = qmethod.getComponentByNumber(i)
base = ms2quant.getIntensityNormalisationBase(gcbn.getName())
component_bases.append(base)
print("Intensity normalisation constants: %s" %(component_bases))
print("Ration normalisation constants: %s" %(rationames))
else:
print("Normalisation: none")
if qmethod.haveOutliers():
q = qmethod.getOutliers()
print("Outliers = %s" %(q.getMethod()))
else:
print("Outliers: none")
print("\n")
proteins = []
for i in ms_range(1, pepsum.getNumberOfHits()):
hit = pepsum.getHit(i)
proteins.append(hit)
j = 0
protein = pepsum.getNextFamilyProtein(i, j)
while protein !=None:
proteins.append(protein)
j += 1
for protein in proteins:
print("%d::%s:" %(protein.getDB(), protein.getAccession()))
for rationame in rationames:
ratio = customquant.getProteinRatio(protein.getAccession(), protein.getDB(), rationame)
if ratio.isMissing():
print("%10s: ---" %(rationame))
else:
print("%10s = %10g, stderr = %10g, N = %d, normal-p = %6g, hyp-p = %6g" %(rationame, ratio.getValue(), ratio.getStandardError_only(), ratio.getSampleSize(), ratio.getNormalityPvalue_only(), ratio.getHypothesisPvalue_only()))
print(" %10s %s" %('Peptide:', rationame))
for i in ms_range(1, protein.getNumPeptides()):
if protein.getPeptideDuplicate(i) == msparser.ms_protein.DUPE_DuplicateSameQuery:
continue
q = protein.getPeptideQuery(i)
p = protein.getPeptideP(i)
peptide = pepsum.getPeptide(q, p)
if not peptide:
continue
values = ''
for rationame in rationames:
quantKey = msparser.ms_peptide_quant_key(q, p)
ratio = ms2quant.getPeptideRatio(quantKey, rationame)
if ratio.isMissing():
values = "---"
continue
if ratio.isInfinite():
values = "###"
else:
values = ratio.getValue()
print(" %-9s %-10s" %("q%d_p%d" %(q, p), str(values)))
print("\n")