From e9a8bf6208f7b5b02789536f0e52e485c44dc473 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Bego=C3=B1a=20Talavera=20And=C3=BAjar?= <begona.talavera@uni.lu> Date: Fri, 17 Mar 2023 09:55:00 +0100 Subject: [PATCH] Upload New File --- patRoon_codes/Suspect_pos_neg.R | 111 ++++++++++++++++++++++++++++++++ 1 file changed, 111 insertions(+) create mode 100644 patRoon_codes/Suspect_pos_neg.R diff --git a/patRoon_codes/Suspect_pos_neg.R b/patRoon_codes/Suspect_pos_neg.R new file mode 100644 index 0000000..021f03a --- /dev/null +++ b/patRoon_codes/Suspect_pos_neg.R @@ -0,0 +1,111 @@ +################################Suspect screening of CSF with patRoon 2.1.0############################### + +###Load library +library(patRoon) + + +###Set all options +#Set the path of the working directory +setwd("../../../../patRoon_codes/") + +#set the name of the output folder +output <- "report_suspect" + +# Modify this value to point to the location of the MetFrag2.4.5-CL.jar file +options(patRoon.path.MetFragCL = "~/opt/metfrag/MetFrag2.4.5-CL.jar") + +options(patRoon.path.obabel = "~/babel/bin/") + +options(patRoon.path.MetFragPubChemLite = ".../patRoon_files/databases/PubChemLite_exposomics_20220729.csv") + + +#Load the csv with the data of the mzML files (note that here only the RP files are shown,as an example but there are two different analysis.csv files por the HILIC analysis) +#Positive ionization files +anaInfoPos <- read.csv("analysis_RP_POS.csv") + +#Negative ionization files +anaInfoNeg <- read.csv("analysis_RP_NEG.csv") + + +#Select the path of the suspect list +suspFile <- read.csv(".../patRoon_files/suspects/name_of_the_suspect_list.csv") + + +###Extract features +library(xcms) + +#To find positive features +fListPos <- findFeatures(anaInfoPos, "xcms3", param = xcms::CentWaveParam(ppm = 5, peakwidth = c(10, 30), mzdiff = -0.001, mzCenterFun = "wMean", integrate = 1L, snthresh = 10, prefilter = + c(3, 100), noise = 7500, fitgauss = FALSE, firstBaselineCheck = TRUE),verbose = TRUE) +#To find negative features + +fListNeg <- findFeatures(anaInfoNeg, "xcms3", param = xcms::CentWaveParam(ppm = 5, peakwidth = c(10, 30), mzdiff = -0.001, mzCenterFun = "wMean", integrate = 1L, snthresh = 10, prefilter = + c(3, 100), noise = 7500, fitgauss = FALSE, firstBaselineCheck = TRUE),verbose = TRUE) +#To merge pos and neg results +fList <- makeSet(fListPos, fListNeg, adducts = c("[M+H]+", "[M-H]-")) + + +###Group features +library(BiocParallel) +BiocParallel::register(BiocParallel::SerialParam(), default = TRUE) + +fGroups <- groupFeatures(fList, "xcms3") + +library ("MetaClean") +fGroups <- filter(fGroups, absMinIntensity = 100, + relMinReplicateAbundance = NULL, maxReplicateIntRSD = NULL, + blankThreshold = 3, removeBlanks = TRUE, + retentionRange = NULL, mzRange = NULL) + +###Screen suspects +fGroups <- screenSuspects(fGroups, suspFile, rtWindow = 12, mzWindow = 0.005, onlyHits = TRUE) + +###Retrieve MS peak lists +avgPListParams <- getDefAvgPListParams(clusterMzWindow = 0.005) +mslists <- generateMSPeakLists(fGroups, "mzr", maxMSRtWindow = 5, + precursorMzWindow = 0.5, avgFeatParams = + avgPListParams, avgFGroupParams = + avgPListParams) + +mslists <- filter(mslists, withMSMS = TRUE, absMSIntThr = NULL, absMSMSIntThr = NULL, relMSIntThr = NULL, + relMSMSIntThr = NULL, topMSPeaks = NULL, topMSMSPeaks = 25, + deIsotopeMS = FALSE, deIsotopeMSMS = FALSE) + + +###Find compound structure candidates +compounds <- generateCompounds(fGroups, mslists, "metfrag", method = "CL", + dbRelMzDev = 5, fragRelMzDev = 5, + fragAbsMzDev = 0.002, + database = "pubchemlite", + scoreTypes = c("fragScore","metFusionScore","score", "individualMoNAScore", + "AnnoTypeCount","Patent_Count", "DisorderDisease", "PubMed_Count"), + maxCandidatesToStop = 100, timeoutRetries = 20) + +###Annotate suspects +fGroups <- annotateSuspects(fGroups, formulas = NULL, + compounds = compounds, MSPeakLists = mslists) + +###Report results + +#Create the csv tables with groups and peak intensities for the set (positive and negative results combined) +reportCSV(fGroups, path = output, reportFeatures = FALSE, + compounds = compounds, compoundsNormalizeScores = "max") + +#Individual features tables +pos_features <- as.data.table(unset(fGroups, "positive")) # table with all feature group data for positive mode +write.csv(pos_features, "report_suspect/pos_features.csv") +neg_features <- as.data.table(unset(fGroups, "negative")) # table with all feature group data for negative mode +write.csv(neg_features, "report_suspect/neg_features.csv") + + +# Summary of MetFrag Results in a a Single Table +MFsummary <- as.data.table(compounds) +outputSummary <- paste(output, "MFsummary.csv", sep = "/") +write.csv(MFsummary, outputSummary) + +#Plots of every group of features (in the set) +reportPDF(fGroups, path = output, reportFGroups = TRUE, + reportFormulaSpectra = TRUE, compounds = compounds, + compoundsNormalizeScores = "max", + MSPeakLists = mslists) + -- GitLab