R/testing.R
In MASHUOA/HiTMaP: A collection of tools for imaging MS data processing
Defines functions
Cardinal_cluster_classification
Pathway_overview_graphite
intensity.colors_customize
Cardinal_old
quiet
Plot_Ion_image_Png2
Plot_Ion_image_Png3
Plot_Ion_image_Png
Kernel_convolution
nice_file_create
fetch_ID_from_DB
Whole_Picture_Cluster_mode
META_pathway_Whole_Picture_discovery_mode
META_pathway_Whole_Picture
META_pathway_ion_image_cardinal
META_pathway_ion_image
META_pathway_search
WorkingDir
simple_ion_image
Generate_HMDBIDS
cluster_image_cardinal
intensity.colors_customize_mono
metwork_cardinal
metwork_maldiquant
Cardinal_utilities
Generatespectrum_workflow_maldiquant
Cardinal
run_spatial_quant_workflow
imaging_default_parameter
RUN_imagine_workflow
# code for testing purpose
RUN_imagine_workflow<-function(){
#Spatial_Quant(ppm=15,adducts=c("M-H","M+Cl"),Quant_list="lipid candidates.csv",cal.mz=T)
#Spatial_Quant(ppm=15,adducts=c("M-H","M+Cl"),Quant_list="lipid candidates.csv",cal.mz=T,Spectrum_feature_summary = T,Region_feature_summary = T,PMF_analysis = F)
imaging_identification(
#==============Choose the imzml raw data file(s) to process
datafile=tk_choose.files(filter = matrix(c( "imzml file", ".imzML",
"Text", ".txt", "All files", "*"),
3, 2, byrow = TRUE),caption = "Choose single or multiple file(s) for analysis"),
threshold=0.05,
ppm=20,
Digestion_site="[G]",
missedCleavages=0:2,
Fastadatabase="murine_matrisome.fasta",
adducts=c("M+H","M+NH4","M+Na"),
#==============TRUE if you want to perform PMF analysis
PMF_analysis=T,
#==============TRUE if you want to generate peptide summary in the Summary folder
Protein_feature_summary=T,
Peptide_feature_summary=T,
#==============TRUE if you want to plot peptide in the Ion images folder, make sure there's imzml file in the folder
plot_ion_image=F,
#==============Set a number if you want a parallel processing
parallel=8,
#==============Set a number (1 to maximum pixels in the data file) if you want to dig more peaks in the raw data
spectra_segments_per_file=5,
#==============True if you want the mean spectrum generated based on spatial spectrum profile
spatialKMeans=T,
Smooth_range=1
)
}
imaging_default_parameter<-function(){
datafile=tk_choose.files(filter = matrix(c( "imzml file", ".imzML",
"Text", ".txt", "All files", "*"),
3, 2, byrow = TRUE),caption = "Choose single or multiple file(s) for analysis")
threshold=0.05
ppm=10
Digestion_site="[G]"
missedCleavages=0:2
Fastadatabase="murine_matrisome.fasta"
adducts=c("M+H","M+NH4","M+Na")
#==============TRUE if you want to perform PMF analysis
PMF_analysis=T
#==============TRUE if you want to generate peptide summary in the Summary folder
Protein_feature_summary=T
Peptide_feature_summary=T
#==============TRUE if you want to plot peptide in the Ion images folder, make sure there's imzml file in the folder
plot_ion_image=T
#==============Set a number if you want a parallel processing
parallel=8
#==============Set a number (1 to maximum pixels in the data file) if you want to dig more peaks in the raw data
spectra_segments_per_file=5
spatialKMeans=T
plot_cluster_image=T
}
run_spatial_quant_workflow<-function(){
imagine_Spatial_Quant(datafile=tk_choose.files(filter = matrix(c( "imzml file", ".imzML",
"Text", ".txt", "All files", "*"),
3, 2, byrow = TRUE),caption = "Choose single or multiple file(s) for analysis"),
threshold=0.00,
ppm=5,
#Fastadatabase="murine_matrisome.fasta",
Quant_list="Metabolites of Interest.csv",
adducts=c("M-H","M+Cl"),
cal.mz=F,
mzlist_bypass=T,
#==============TRUE if you want to plot protein PMF result
PMF_analysis=TRUE,
#==============TRUE if you want to generate peptide summary in the Summary folder
Protein_feature_summary=F,
Peptide_feature_summary=F,
#PLOT_PMF_Protein=FALSE,
#==============TRUE if you want to plot peptide in the Ion images folder, make sure there's imzml file in the folder
plot_ion_image=FALSE,
#==============Set a number if you want a parallel processing
parallel=detectCores()/2,
#==============Set a number (1 to maximum pixels in the data file) if you want to dig more peaks in the raw data
spectra_segments_per_file=4,
spatialKMeans=F,
Smooth_range=1,
Virtual_segmentation=T,
Virtual_segmentation_rankfile="Z:\\George skyline results\\maldiimaging\\Maldi_imaging - Copy\\radius_rank_bovin.csv",
Spectrum_feature_summary=F,
Region_feature_summary=F,
Region_feature_analysis=F,
Cluster_level=="Low")
#workflow for aged human lens
imagine_Spatial_Quant(
#==============Choose the imzml raw data file(s) to process make sure the fasta file in the same folder
datafile=tk_choose.files(filter = matrix(c( "imzml file", ".imzML",
"Text", ".txt", "All files", "*"),
3, 2, byrow = TRUE),caption = "Choose single or multiple file(s) for analysis"),
threshold=0.00,
ppm=20,
#Fastadatabase="murine_matrisome.fasta",
Quant_list="lipid candidates.csv",
adducts=c("M-H","M+Cl"),
cal.mz=T,
mzlist_bypass=F,
#==============TRUE if you want to plot protein PMF result
PMF_analysis=TRUE,
#==============TRUE if you want to generate protein summary in the Summary folder
Protein_feature_summary=T,
#==============TRUE if you want to generate protein cluster image in the Summary folder
plot_cluster_image=T,
Peptide_feature_summary=T,
#PLOT_PMF_Protein=FALSE,
#==============TRUE if you want to plot peptide in the Ion images folder, make sure there's imzml file in the folder
plot_ion_image=FALSE,
#==============Set a number if you want a parallel processing
parallel=detectCores(),
#==============Set a number (1 to maximum pixels in the data file) if you want to dig more peaks in the raw data
spectra_segments_per_file=5,
spatialKMeans=F,
Smooth_range=1,
Virtual_segmentation=T,
Virtual_segmentation_rankfile="Z:\\George skyline results\\maldiimaging\\Maldi_imaging - Copy\\radius_rank.csv",
Spectrum_feature_summary=T,
Region_feature_summary=T,
Region_feature_analysis=T,
plot_each_metabolites=F,
Cluster_level="High",
Region_feature_analysis_bar_plot=F,
norm_datafiles=T,
norm_Type="Median"
)
imagine_Spatial_Quant(
#==============Choose the imzml raw data file(s) to process make sure the fasta file in the same folder
datafile=tk_choose.files(filter = matrix(c( "imzml file", ".imzML",
"Text", ".txt", "All files", "*"),
3, 2, byrow = TRUE),caption = "Choose single or multiple file(s) for analysis"),
threshold=0.00,
ppm=25,
#Fastadatabase="murine_matrisome.fasta",
Quant_list="lipid blast candidates.csv",
adducts=c("M-H","M+Cl"),
cal.mz=F,
mzlist_bypass=F,
#==============TRUE if you want to plot protein PMF result
PMF_analysis=TRUE,
#==============TRUE if you want to generate protein summary in the Summary folder
Protein_feature_summary=T,
#==============TRUE if you want to generate protein cluster image in the Summary folder
plot_cluster_image=T,
Peptide_feature_summary=T,
#PLOT_PMF_Protein=FALSE,
#==============TRUE if you want to plot peptide in the Ion images folder, make sure there's imzml file in the folder
plot_ion_image=FALSE,
#==============Set a number if you want a parallel processing
parallel=detectCores(),
#==============Set a number (1 to maximum pixels in the data file) if you want to dig more peaks in the raw data
spectra_segments_per_file=5,
spatialKMeans=F,
Smooth_range=1,
Virtual_segmentation=T,
Virtual_segmentation_rankfile="Z:\\George skyline results\\maldiimaging\\Maldi_imaging - Copy\\radius_rank.csv",
Spectrum_feature_summary=T,
Region_feature_summary=T,
Region_feature_analysis=T,
plot_each_metabolites=F,
Cluster_level="High",
Region_feature_analysis_bar_plot=F,
norm_datafiles=T
)
}
Cardinal<- function(){
workdir <- WorkingDir()
datafile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.imzML$")
listfile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.csv$")
datafile<-gsub(".imzML", "", datafile)
name <-gsub(paste(workdir,"\\/",sep=""),"",datafile[3])
#imdata <- readImzML(name, folder, as="MSImagingExperiment")
imdata <- readImzML(name, folder, attach.only=TRUE, as="MSImagingExperiment")
summarize(data, sum, .by="pixel")
tmp <- imdata %>%
smoothSignal() %>%
reduceBaseline() %>%
peakPick() %>%
peakFilter() %>%
select(x == 1, y == 1)
process(plot=TRUE,
par=list(layout=c(1,3)),
BPPARAM=SerialParam())
plot(data, pixel=1)
plot(data, coord=list(x=2, y=2))
pattern <- factor(c(0, 0, 2, 2, 0, 0, 0, 0, 0, 0, 2, 2, 0,
0, 0, 0, 0, 0, 0, 1, 2, 2, 0, 0, 0, 0, 0, 2, 1, 1, 2,
2, 0, 0, 0, 0, 0, 1, 2, 2, 2, 2, 0, 0, 0, 0, 1, 2, 2,
2, 2, 2, 0, 0, 0, 0, 2, 2, 2, 2, 2, 2, 2, 0, 0, 0, 2,
2, 0, 0, 0, 0, 0, 0, 2, 2, 0, 0, 0, 0, 0),
levels=c(0,1,2), labels=c("blue", "black", "red"))
tmp <- imdata %>%
smoothSignal() %>%
reduceBaseline() %>%
peakPick() %>%
peakFilter() %>%
select(x == 1, y == 1) %>%
process(plot=TRUE,
par=list(layout=c(1,3)),
BPPARAM=SerialParam())
}
Generatespectrum_workflow_maldiquant<-function(){
table(sapply(MALDI_IMAGE, length))
plot(MALDI_IMAGE[[1]])
workdir <- WorkingDir()
datafile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.imzML$")
datafile<-gsub(".imzML", "", datafile)
MALDI_IMAGE <- quiet(importImzMl(paste(file.path( datafile[i]),".imzml",sep="")))
#use the square root transformation to simplify graphical visualization and to overcome the potential dependency of the variance from the mean.
spectra <- transformIntensity(MALDI_IMAGE, method="sqrt")
#use a 21 point Savitzky-Golay-Filter (Savitzky and Golay, 1964) to smooth the spectra.
spectra <- smoothIntensity(spectra, method="SavitzkyGolay", halfWindowSize=10)
#Before we correct the baseline we visualize it. Here we use the SNIP algorithm (Ryan et al., 1988).
baseline <- estimateBaseline(spectra[[16]], method="SNIP", iterations=100)
plot(spectra[[16]])
lines(baseline, col="red", lwd=2)
#If we are satisfied with our estimated baseline we remove it.
spectra <- removeBaseline(spectra, method="SNIP",iterations=100)
plot(spectra[[1]])
#spectra <- calibrateIntensity(spectra, method="TIC")
#do the alignment if there's a huge mass shift between each spectrum
spectra <- alignSpectra(spectra,halfWindowSize=50,SNR=3,tolerance=0.05,warpingMethod="lowess")
samples <- factor(sapply(spectra, function(x)metaData(x)$sampleName))
avgSpectra <- averageMassSpectra(spectra, labels=samples, method="mean")
noise <- estimateNoise(avgSpectra[[1]])
plot(avgSpectra[[1]], xlim=c(500, 4000), ylim=c(0, 0.002))
lines(noise, col="red")
lines(noise[,1], noise[, 2]*2, col="blue")
peaks <- detectPeaks(avgSpectra, method="MAD",halfWindowSize=20, SNR=2)
plot(avgSpectra[[1]], xlim=c(4000, 5000), ylim=c(0, 0.002))
points(peaks[[1]], col="red", pch=4)
peaks <- binPeaks(peaks, tolerance=0.002)
peaks <- filterPeaks(peaks, minFrequency=0.25)
}
Cardinal_utilities<- function(){
folder <- WorkingDir()
datafile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.imzML$")
listfile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.csv$")
datafile<-gsub(".imzML", "", datafile)
name <-gsub(paste(folder,"\\/",sep=""),"",datafile[1])
#imdata <- readImzML(name, folder, as="MSImagingExperiment")
imdata <- readImzML(name, folder, attach.only=TRUE, as="MSImagingExperiment")
summarize(data, sum, .by="pixel")
tmp <- data %>%
smoothSignal() %>%
reduceBaseline() %>%
peakPick() %>%
peakFilter() %>%
select(x == 1, y == 1)
process(plot=TRUE,
par=list(layout=c(1,3)),
BPPARAM=SerialParam())
plot(data, pixel=1)
plot(data, coord=list(x=2, y=2))
pattern <- factor(c(0, 0, 2, 2, 0, 0, 0, 0, 0, 0, 2, 2, 0,
0, 0, 0, 0, 0, 0, 1, 2, 2, 0, 0, 0, 0, 0, 2, 1, 1, 2,
2, 0, 0, 0, 0, 0, 1, 2, 2, 2, 2, 0, 0, 0, 0, 1, 2, 2,
2, 2, 2, 0, 0, 0, 0, 2, 2, 2, 2, 2, 2, 2, 0, 0, 0, 2,
2, 0, 0, 0, 0, 0, 0, 2, 2, 0, 0, 0, 0, 0),
levels=c(0,1,2), labels=c("blue", "black", "red"))
tmp <- imdata %>%
smoothSignal() %>%
reduceBaseline() %>%
peakPick() %>%
peakFilter() %>%
select(x == 1, y == 1) %>%
process(plot=TRUE,
par=list(layout=c(1,3)),
BPPARAM=SerialParam())
}
####load pathway database
#Map_SMP_MP <- read.csv("DB//smpdb_pathways.csv",as.is = TRUE)
####load raw data annotation
####load adducts data
#adducts<-read.csv("DB//adducts.csv",as.is = TRUE)
####load SMPDB data
#SMPDB <- read.csv(file="DB/SMPdb.csv",stringAsFactors=FALSE,as.is = TRUE)
#SMPDBFREQ <- read.csv(file="DB/SMPdb_ID_FREQ.csv",col.names= c("","SMPDB.ID","Freq"),as.is = TRUE)
##work dir
##load raw data file
#MALDI_IMAGE <- importImzMl(file.path( "NEDC.imzML"))
######Example steps
metwork_maldiquant<- function(workdir= tk_choose.dir(caption = "Select working directory")){
datafile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.imzML$")
datafile<-gsub(".imzML", "", datafile)
setwd("Z:\\George skyline results\\maldiimaging\\Maldi_imaging - Copy")
HMDBids <- read.csv("combine.csv", fill = TRUE,skip=2,as.is = TRUE)
SMPDB <- read.csv(file="DB/SMPdb.csv",as.is = TRUE)
for (i in 1: length(datafile)){
if (dir.exists(datafile[i])==FALSE){dir.create(datafile[i])}
HMDBDATA <- Generate_HMDBIDS(HMDBids)
SMPLIST <- META_pathway_search(HMDBDATA, SMPDB, datafile[i], PWtype="Metabolic")
SMPRESULT <-read.csv(file=paste(datafile[i],sep="","/","SMPRESULT.csv"),as.is = TRUE)
options(warn=-1)
MALDI_IMAGE <- quiet(importImzMl(paste(file.path( datafile[i]),".imzml",sep="")))
options(warn=0)
META_pathway_ion_image(SMPLIST, SMPRESULT, datafile[i] ,MALDI_IMAGE, Image_Type="",PWpScore=0.26,plotTolerance=0.125 )
META_pathway_Whole_Picture(SMPLIST, SMPRESULT, datafile[i], PWpScore=0.26,pathway_mode="")
META_pathway_Whole_Picture_discovery_mode(SMPLIST, SMPRESULT, datafile[i], SMPDB, adductslist=c("M-H","M+Cl"), PWpScore=0.26)
}
}
metwork_cardinal<- function(workdir= tk_choose.dir(caption = "Select working directory")){
datafile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.imzML$")
datafile<-gsub(".imzML", "", datafile)
setwd("Z:\\George skyline results\\maldiimaging\\Maldi_imaging - Copy")
HMDBids <- read.csv("combine.csv", fill = TRUE,skip=2,as.is = TRUE)
SMPDB <- read.csv(file="DB/SMPdb.csv",as.is = TRUE)
for (i in 1: length(datafile)){
if (dir.exists(datafile[i])==FALSE){dir.create(datafile[i])}
HMDBDATA <- Generate_HMDBIDS(HMDBids)
SMPLIST <- META_pathway_search(HMDBDATA, SMPDB, datafile[i], PWtype="Metabolic")
SMPRESULT <-read.csv(file=paste(datafile[i],sep="","/","SMPRESULT.csv"),as.is = TRUE)
options(warn=-1)
imdata <- quiet(Load_Cardinal_imaging(paste(file.path( datafile[i]),sep=""),resolution = 20))
options(warn=0)
META_pathway_ion_image_cardinal(SMPLIST, SMPRESULT, datafile[i] ,imdata, Image_Type="",PWpScore=0.26,plotTolerance=20 ,)
filename=gsub(dirname(datafile[i]),"",datafile[i])
for (SMPDB.ID in SMPRESULT$SMPDB.ID[SMPRESULT$pScore>=PWpScore]){
cluster_image_cardinal(SMPDB.ID,imdata,SMPLIST,SMPRESULT$Name[SMPRESULT$SMPDB.ID==SMPDB.ID],ppm=ppm)
#cluster_image_cardinal(SMPDB.ID,imdata_istd,SMPLIST,SMPRESULT$Name[SMPRESULT$SMPDB.ID==SMPDB.ID])
}
#lapply(SMPRESULT$SMPDB.ID[SMPRESULT$pScore>=PWpScore],cluster_image_cardinal,imdata,SMPLIST)
META_pathway_Whole_Picture(SMPLIST, SMPRESULT, datafile[i], PWpScore=0.26,pathway_mode="")
META_pathway_Whole_Picture_discovery_mode(SMPLIST, SMPRESULT, datafile[i], SMPDB, adductslist=c("M-H","M+Cl"), PWpScore=0.26)
}
}
##get HMDB IDs list for pathway search
intensity.colors_customize_mono <- function(n = 100, alpha = 1,colorstart="transparent",colorend="red") {
col1 <- colorRamp(colors=c(colorstart, colorend),interpolate = "linear")
col1
}
cluster_image_cardinal<-function(clusterID,imdata,SMPLIST,clustername,ppm=20){
library(RColorBrewer)
outputpng=paste(getwd(),"/",clustername,"_cluster_plot",'.png',sep="")
candidate=SMPLIST[SMPLIST$SMPDB.ID==clusterID,]
candidateunique=unique(candidate[,"mz"])
mycol <- RColorBrewer::brewer.pal(length(candidateunique),"Set1")
png(outputpng,width = 720,height = 720, bg = "transparent")
image(imdata, mz=candidateunique,
col=mycol,
#contrast.enhance = contrast.enhance,
smooth.image = smooth.image ,
superpose=TRUE,normalize.image="linear")
dev.off()
pngfile<-image_read(outputpng)
pngfile<-image_border(pngfile, "transparent", "30x30")
pngfile<-image_annotate(pngfile,paste(clustername),gravity = "south",size = 40,color = "white")
pngfile<-image_trim(pngfile)
image_write(pngfile,outputpng)
for (i in 1:length(candidateunique)){
outputpngsum=paste(getwd(),"/",clustername,"_cluster_plot_",candidateunique[i],'.png',sep="")
png(outputpngsum,width = 720,height = 720, bg = "transparent")
par(oma=c(0, 0, 0, 0),tcl = NA,mar=c(0, 0, 1, 1),mfrow = c(1, 1),
bty="n",pty="s",xaxt="n",
yaxt="n",
no.readonly = TRUE,ann=FALSE)
#image(imdata, mz=candidateunique[i], col=mycol[i], superpose=F,normalize.image="linear")
col.regions <- gradient.colors(100, start="transparent", end=mycol[i])
image(imdata, mz=candidateunique[i],
contrast.enhance=contrast.enhance,
smooth.image = smooth.image ,
col.regions=col.regions,normalize.image="linear")
dev.off()
}
}
#plot cluster image all in one file
Generate_HMDBIDS <- function(HMDBids){
HMDBDATA<-NULL
HMDBSearchlist<-str_split(paste(HMDBids$moleculeIds,collapse=", "),pattern=", ")[[1]]
Formular<-str_split(paste(HMDBids$moleculeNames,collapse=", "),pattern=", ")[[1]]
HMDBDATA$HMDB.ID<-HMDBSearchlist
HMDBDATA$moleculeNames<-Formular
MZLIST <- NULL
adduct <- NULL
options(warn=-1)
for (i in 1: length(HMDBSearchlist)){
MZLIST[i]<-as.numeric(HMDBids[grep(HMDBSearchlist[i],HMDBids[,"moleculeIds"]),"mz"])
adduct[i]<-as.character(HMDBids[grep(HMDBSearchlist[i],HMDBids[,"moleculeIds"]),"adduct"])
}
HMDBDATA$mz<-MZLIST
HMDBDATA$adduct<-adduct
options(warn=0)
HMDBDATA<-as.data.frame(HMDBDATA)
HMDBDATA<-unique(HMDBDATA)
return(HMDBDATA)
}
##simple_ion_image: define the work dir read all csv and enerate the ion image respectively
simple_ion_image<-function(workdir=WorkingDir(),
Image_Type="",
plotTolerance=5 ,
creat_new_file=TRUE,
Denoising=0,
color = "black",
interpolate =FALSE,
...){
datafile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.imzML$")
listfile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.csv$")
masslist<-NULL
for (z in 1:length(datafile)){
datafile[z]<-gsub(".imzML", "", datafile[z])
options(warn=-1)
MALDI_IMAGE <- importImzMl(paste(file.path(datafile[z]),".imzML",sep=""))
options(warn=0)
if (dir.exists(datafile[z])==FALSE){dir.create(datafile[z])}
for (i in 1:length(listfile)){
foldername<-str_remove(listfile[i],paste(workdir,"/",sep=""))
foldername<-str_remove(foldername,".csv")
outputDir<-paste(datafile[z],"/",foldername,sep="")
if (dir.exists(outputDir)==FALSE){dir.create(outputDir)}
masslist<-read.csv(listfile[i],header = TRUE)
for (j in 1:length(masslist$moleculeNames)){
if (is.na(masslist$mz[j])==FALSE){
for (k in 1:length(str_split(masslist$moleculeNames[j],", ")[[1]])){
Plot_Ion_image_Png3(outputDir, MALDI_IMAGE, windows_filename(str_split(masslist$moleculeNames[j],", ")[[1]][k]), masslist$mz[j],masslist$adduct[j], Tolerance= plotTolerance*masslist$mz[j]/1000000, title="",Creat_new_file=TRUE,Neighbour = Denoising,color = color,interpolate =interpolate)
}}}
}
}
}
WorkingDir <-function(Mainfolder = tk_choose.dir(caption = "Select working directory")){
Mainfolder <- Mainfolder
setwd(Mainfolder)
Mainfolder
}
##Pathway search return search result
META_pathway_search <- function(HMDBDATA, SMPDB, WKdir, PWtype="Metabolic", ...) {
SMPLIST<- NULL
SMPMZlist<-NULL
for (i in 1: length(HMDBDATA[,"HMDB.ID"])){
SMPquery<-NULL
SMPquery<-SMPDB[`&`(SMPDB[,"Pathway.Subject"]==PWtype,SMPDB[,"HMDB.ID"] == HMDBDATA[,"HMDB.ID"][i]),"SMPDB.ID"]
SMPLIST$SMPDB.ID<-c(SMPLIST$SMPDB.ID,SMPquery)
SMPLIST$HMDB.ID<-c(SMPLIST$HMDB.ID, rep(as.character(HMDBDATA[,"HMDB.ID"][i]),length(SMPquery)))
SMPLIST$mz<-c(SMPLIST$mz, rep(HMDBDATA[,"mz"][i],length(SMPquery)))
SMPLIST$moleculeNames<-c(SMPLIST$moleculeNames, rep(as.character(HMDBDATA[,"moleculeNames"][i]),length(SMPquery)))
SMPLIST$KEGG.ID<-c(SMPLIST$KEGG.ID, rep(SMPDB[SMPDB[,"HMDB.ID"]==HMDBDATA[,"HMDB.ID"][i],"KEGG.ID"][1],length(SMPquery)))
SMPLIST$Metabolite.ID<-c(SMPLIST$Metabolite.ID, rep(SMPDB[SMPDB[,"HMDB.ID"]==HMDBDATA[,"HMDB.ID"][i],"Metabolite.ID"][1],length(SMPquery)))
SMPLIST$adduct<-c(SMPLIST$adduct, rep(as.character(HMDBDATA[HMDBDATA[,"HMDB.ID"]==HMDBDATA[,"HMDB.ID"][i],"adduct"][1]),length(SMPquery)))
}
if (dir.exists(paste(WKdir,sep="","/"))==FALSE){dir.create(paste(WKdir,sep="","/"))}
SMPLIST<-as.data.frame(SMPLIST)
write.csv(SMPLIST,paste(WKdir,sep="","/","SMPqueryRESULT.csv"),row.names=FALSE)
SMPRESULT<-as.data.frame(table(SMPLIST$SMPDB.ID))
write.csv(SMPRESULT,paste(WKdir,sep="","/","SMPRESULT.csv"),row.names=FALSE)
SMPRESULT<-read.csv(file=paste(WKdir,sep="","/","SMPRESULT.csv"),col.names= c("SMPDB.ID","Freq"))
SMPRESULTfreq<-NULL
for (i in 1: length(SMPRESULT$SMPDB.ID)){
SMPRESULTfreq<-c(SMPRESULTfreq,SMPDBFREQ[SMPDBFREQ[,"SMPDB.ID"]==as.character(SMPRESULT[i,"SMPDB.ID"]),"Freq"])
SMPRESULT$Name[i]<-as.character(Map_SMP_MP[Map_SMP_MP[,"SMPDB.ID"] %in% SMPRESULT[,"SMPDB.ID"][i],"Name"])
}
SMPRESULT$DBfreq<-SMPRESULTfreq
SMPRESULT$pScore<-SMPRESULT$Freq/SMPRESULT$DBfreq
write.csv(SMPRESULT,paste(WKdir,sep="","/","SMPRESULT.csv"),row.names=FALSE)
return(SMPLIST)
}
###Ion Image manipulation
META_pathway_ion_image<- function(SMPLIST, SMPRESULT, datafile,MALDI_IMAGE, Image_Type="",PWpScore=0.26,plotTolerance=0.125 ,creat_new_file=TRUE,...) {
if (dir.exists(datafile)==FALSE){dir.create(datafile)}
for (j in 1:length(SMPRESULT$SMPDB.ID)){
if (SMPRESULT$pScore[j]>=PWpScore){
WKdir<-paste(datafile,"/",SMPRESULT$SMPDB.ID[j],sep="")
PWSUBLIST<-SMPLIST[SMPLIST$SMPDB.ID==SMPRESULT$SMPDB.ID[j],]
if (dir.exists(WKdir)==FALSE){dir.create(WKdir)}
for (i in 1:length(PWSUBLIST$mz)){
Plot_Ion_image_Png(WKdir,MALDI_IMAGE,PWSUBLIST$moleculeNames[i],PWSUBLIST$mz[i],PWSUBLIST$adduct[i],Tolerance= plotTolerance, title=Image_Type,Creat_new_file=creat_new_file)
}}
}}
META_pathway_ion_image_cardinal<- function(SMPLIST, SMPRESULT, datafile,MALDI_IMAGE, Image_Type="",PWpScore=0.26,plotTolerance=0.125 ,creat_new_file=TRUE,...) {
if (dir.exists(datafile)==FALSE){dir.create(datafile)}
for (j in 1:length(SMPRESULT$SMPDB.ID)){
if (SMPRESULT$pScore[j]>=PWpScore){
WKdir<-paste(datafile,"/",SMPRESULT$SMPDB.ID[j],sep="")
PWSUBLIST<-SMPLIST[SMPLIST$SMPDB.ID==SMPRESULT$SMPDB.ID[j],]
if (dir.exists(WKdir)==FALSE){dir.create(WKdir)}
for (i in 1:length(PWSUBLIST$mz)){
Plot_Ion_image_Png_Cardinal(WKdir,MALDI_IMAGE,PWSUBLIST$moleculeNames[i],PWSUBLIST$mz[i],PWSUBLIST$adduct[i],Tolerance= plotTolerance, title=Image_Type,Creat_new_file=creat_new_file,contrast.enhance = "suppression")
}}
}}
###DRAW META_pathway_Whole_Picture
META_pathway_Whole_Picture<- function(SMPLIST, SMPRESULT, datafile, PWpScore=0.26,pathway_mode="", ...) {
MPID <-NULL
if (dir.exists(datafile)==FALSE){dir.create(datafile)}
for (j in 1:length(SMPRESULT$SMPDB.ID)){
if (SMPRESULT$pScore[j]>=PWpScore){
WKdir<-paste(datafile,"/",SMPRESULT$SMPDB.ID[j],sep="")
if (exists("Map_SMP_MP")==FALSE){Map_SMP_MP <- read.csv("DB/smpdb_pathways.csv")}
MPID <- sprintf("PW%06d", match(SMPRESULT$SMPDB.ID[j],Map_SMP_MP[,1]))
#PWsvg <- htmlParse(paste("DB/smpdb_svg/",MPID ,".svg",sep=""))
PWSUBLIST<-SMPLIST[SMPLIST$SMPDB.ID==SMPRESULT$SMPDB.ID[j],c("HMDB.ID","mz","moleculeNames","Metabolite.ID")]
PW_svg_line<- readLines(paste("DB/smpdb_svg/",MPID ,".svg",sep=""))
for (i in 1:length(PWSUBLIST$mz)){
pngfile<- NULL
res <- try(pngfile<-image_read(paste(WKdir,"/",PWSUBLIST$moleculeNames[i],'.png',sep="")),silent = TRUE)
if (class(res) != "try-error"){
Ion_image_svg_base64 <- base64Encode(readBin(paste(WKdir,"/",PWSUBLIST$moleculeNames[i],'.png',sep=""), "raw", file.info(paste(WKdir,"/",PWSUBLIST$moleculeNames[i],'.png',sep=""))[1, "size"]), "txt")
res<-try(Start_moldb_cache<-grep(paste('id=','"',PWSUBLIST$Metabolite.ID[i],'_moldb_cache',sep=""),as.character(PW_svg_line),value=FALSE),silent = TRUE)
res2<-try(End_moldb_cache<-grep("pointer-events=",as.character(PW_svg_line[Start_moldb_cache:length(PW_svg_line)]),value=FALSE)[1],silent = TRUE)
if ("&"(class(res) != "try-error",class(res2) != "try-error")){
PW_svg_line<-c(PW_svg_line[1:Start_moldb_cache], paste("<image xlink:href=\"data:image/png;base64,",Ion_image_svg_base64,aep=""),PW_svg_line[(Start_moldb_cache+End_moldb_cache-1):length(PW_svg_line)])
if (pathway_mode=="D-mode"){write(PW_svg_line,file=paste(datafile,"/",SMPRESULT$SMPDB.ID[j]," D-mode.svg",sep=""))}else
{write(PW_svg_line,file=paste(datafile,"/",SMPRESULT$SMPDB.ID[j],".svg",sep=""))}}
}
}
}
}
}
#Metapathway Image discovery_mode_list
META_pathway_Whole_Picture_discovery_mode<- function(SMPLIST, SMPRESULT, datafile, SMPDB, adductslist="M-H", PWpScore=0.26, ...) {
MPID <-NULL
SMPLIST_D_MODE<-SMPDB[`&`(SMPDB[,"SMPDB.ID"] %in% SMPRESULT[SMPRESULT[,"pScore"]>0.26,"SMPDB.ID"],!(SMPDB[,"HMDB.ID"] %in% unique(SMPLIST[,"HMDB.ID"]))),]
SMPLIST_D_MODE$mass<- t(as.data.frame(lapply(SMPLIST_D_MODE$Formula, getMonomass)))
SMPLIST_D_MODE$moleculeNames<-SMPLIST_D_MODE$Metabolite.Name
for (adduction in 1:length(adductslist)){
SMPLIST_D_MODE$adduct<-adductslist[adduction]
SMPLIST_D_MODE$mz <- SMPLIST_D_MODE$mass+ adducts[adducts[,"Name"]==adductslist[adduction],"Mass"]
if (adduction==1){
META_pathway_ion_image(SMPLIST_D_MODE, SMPRESULT, datafile,MALDI_IMAGE, Image_Type="D-mode", PWpScore=0.26,plotTolerance=0.125,creat_new_file=TRUE)
}else{
META_pathway_ion_image(SMPLIST_D_MODE, SMPRESULT, datafile,MALDI_IMAGE, Image_Type="D-mode", PWpScore=0.26,plotTolerance=0.125,creat_new_file=FALSE)
}
}
options(warn=-1)
SMPLIST_D_MODE_combine<-bind_rows(SMPLIST_D_MODE[,colnames(SMPLIST)],SMPLIST)
options(warn=0)
META_pathway_Whole_Picture(SMPLIST_D_MODE_combine, SMPRESULT, datafile, PWpScore=0.26,pathway_mode="D-mode")
}
#Cluster_mode_list
Whole_Picture_Cluster_mode<- function(workdir=WorkingDir(),
cluster_list_file="Cluster.csv",
clustercol="Protein",
Sample_list_file="SMAPLE.csv",
adductslist="M-H",
Threshold=0.15,...
) {
cluster_list<-read.csv(cluster_list_file)
Sample_list<-read.csv(Sample_list)
for (i in 1:length(Sample_list$Sample)){
MALDI_IMAGE_list[as.numeric(Sample_list$Order[i])] <- importImzMl(paste(file.path(Sample_list$Sample[i]),".imzml",sep=""))
}
MPID <-NULL
SMPLIST_D_MODE<-SMPDB[`&`(SMPDB[,"SMPDB.ID"] %in% SMPRESULT[SMPRESULT[,"pScore"]>0.26,"SMPDB.ID"],!(SMPDB[,"HMDB.ID"] %in% unique(SMPLIST[,"HMDB.ID"]))),]
SMPLIST_D_MODE$mass<- t(as.data.frame(lapply(SMPLIST_D_MODE$Formula, getMonomass)))
SMPLIST_D_MODE$moleculeNames<-SMPLIST_D_MODE$Metabolite.Name
for (adduction in 1:length(adductslist)){
SMPLIST_D_MODE$adduct<-adductslist[adduction]
SMPLIST_D_MODE$mz <- SMPLIST_D_MODE$mass+ adducts[adducts[,"Name"]==adductslist[adduction],"Mass"]
if (adduction==1){
META_pathway_ion_image(SMPLIST_D_MODE, SMPRESULT, datafile,MALDI_IMAGE, Image_Type="D-mode", PWpScore=0.26,plotTolerance=0.125,creat_new_file=TRUE)
}else{
META_pathway_ion_image(SMPLIST_D_MODE, SMPRESULT, datafile,MALDI_IMAGE, Image_Type="D-mode", PWpScore=0.26,plotTolerance=0.125,creat_new_file=FALSE)
}
}
options(warn=-1)
SMPLIST_D_MODE_combine<-bind_rows(SMPLIST_D_MODE[,colnames(SMPLIST)],SMPLIST)
options(warn=0)
META_pathway_Whole_Picture(SMPLIST_D_MODE_combine, SMPRESULT, datafile, PWpScore=0.26,pathway_mode="D-mode")
}
#for (j in 1:length(SMPRESULT$SMPDB.ID)){
#if (SMPRESULT$pScore[j]>=PWpScore){
#dir<-paste(datafile,"/",SMPRESULT$SMPDB.ID[j],sep="")
#discovery_moldb_list<-NULL
#if (exists("Map_SMP_MP")==FALSE){Map_SMP_MP <- read.csv("DB/smpdb_pathways.csv")}
#MPID <- sprintf("PW%06d", match(SMPRESULT$SMPDB.ID[j],Map_SMP_MP[,1]))
#PWsvg <- htmlParse(paste("DB/smpdb_svg/",MPID ,".svg",sep=""))
#PWSUBLIST<-SMPLIST[SMPLIST$SMPDB.ID==SMPRESULT$SMPDB.ID[j],c("HMDB.ID","mz","moleculeNames","Metabolite.ID")]
#PW_svg_line<- readLines(paste("DB/smpdb_svg/",MPID ,".svg",sep=""))
#discovery_moldb_cache<-grep('<g id="location" transform=',as.character(PW_svg_line),value=FALSE)
#discovery_moldb_cache<-str_locate(PW_svg_line[discovery_moldb_cache],"data-element-id=")
#for (i in 1 : length(discovery_moldb_cache)){
#temp<-substr(PW_svg_line[discovery_moldb_cache[i]],discovery_moldb_startpos[i,"end"]+2,discovery_moldb_startpos[i,"end"]+2+9)
#discovery_moldb_cache<-grep('data-element-type="compound" data-element-id=',as.character(PW_svg_line),value=FALSE)
#discovery_moldb_startpos<-str_locate(PW_svg_line[discovery_moldb_cache],"data-element-id=")
#discovery_moldb_list<-substr(PW_svg_line[discovery_moldb_cache],discovery_moldb_startpos[,"end"]+2,discovery_moldb_startpos[,"end"]+2+9)
#discovery_moldb_list <- union(discovery_moldb_list,discovery_moldb_list)
#discovery_moldb_list<- fetch_ID_from_DB(discovery_moldb_list,SMPDB,SMPRESULT$SMPDB.ID[j])
##fetch_ID_from SMPDB
fetch_ID_from_DB <- function(Fetch_id_list, SMPDB, SMPDB.ID){
if (SMPDB.ID != "") {for (z in 1: length(Fetch_id_list)){temp3[z,]<-subset(SMPDB , Metabolite.ID == as.character(Fetch_id_list)[z] & SMPDB.ID ==SMPRESULT$SMPDB.ID[j])}}else{
for (z in 1: length(Fetch_id_list)){temp3[z,]<-subset(SMPDB , Metabolite.ID == as.character(Fetch_id_list)[z] )}}
return(temp3)
}
nice_file_create <- function(filename){
ok <- file.create(filename)
if(!ok)
{return(windows_filename2(filename))
}else{return(filename)}
}
Kernel_convolution<-function(Neighbour,resolution=1){
Neighbour<-1+(Neighbour*2*resolution)
if (Neighbour==1){kern = matrix(1, ncol = Neighbour, nrow = Neighbour)}else{
kern = matrix(1/(Neighbour*Neighbour-1), ncol = Neighbour, nrow = Neighbour)
kern[(Neighbour+1)/2,(Neighbour+1)/2]=0}
kern
}
## Basic Function PNG ion image plot
Plot_Ion_image_Png<- function(WKdir, imagefile, png_filename, mz,adduct="M-H", Tolerance= 0.25, title="",Neighbour = 0,Creat_new_file=TRUE,color="black",interpolate =FALSE){
#x11()
#dev.copy(png,paste(dir,"/",png_filename,'.png',sep=""))
Tolerance<-round(Tolerance,digits = 5)
pngfillewrite<-paste(WKdir,"/",png_filename,'.png',sep="")
#pngfillewrite<-nice_file_create(pngfillewrite)
png(paste(WKdir,"/","temp.png",sep=""), bg = "transparent")
try(plotMsiSlice(imagefile ,mz , tolerance=Tolerance, legend=FALSE,colRamp=colorRamp(c("black", "blue", "green", "yellow", "red","#FF00FF","white")),interpolate =interpolate ),silent = TRUE)
dev.off()
try(pngfile<-image_read(paste(WKdir,"/","temp.png",sep="")))
res <- try(pngfile<-image_trim(pngfile),silent = TRUE)
if (class(res) != "try-error"){
#resolution<-imagefile[[1]]@metaData[["imaging"]][["Size"]]["x"]
kern=Kernel_convolution(Neighbour,resolution)
if (title=="D-mode"){
pngfile<-image_border(pngfile, "transparent", "0x30")
pngfile<-image_annotate(pngfile,paste(png_filename,"D-mode","\n",adduct,sprintf("%g",mz),"±",Tolerance),gravity = "south",size = 9)
pngfile<-image_trim(pngfile)
res2 <- try(pngfileoriginal<-image_read(pngfillewrite),silent = TRUE)
if ('&'(class(res2) != "try-error",Creat_new_file==FALSE)){
pngfile<-image_append(c(pngfileoriginal,pngfile))
image_write(pngfile,pngfillewrite)}
if ('&'(class(res2) != "try-error",Creat_new_file==TRUE)){
image_write(pngfile,pngfillewrite)}
if (class(res2) == "try-error"){
image_write(pngfile,pngfillewrite)}
}else{
pngfile<-image_convolve(pngfile, kern)
pngfile<-image_border(pngfile, "transparent", "20x20")
pngfile<-image_annotate(pngfile,paste(png_filename,"\n", adduct,sprintf("%g",mz),"±",Tolerance),gravity = "south",size = 10,color = color)
pngfile<-image_trim(pngfile)
image_write(pngfile,pngfillewrite)
}
}
file.remove(paste(WKdir,"/","temp.png",sep=""))
}
Plot_Ion_image_Png3<- function(WKdir, imagefile, png_filename, mz,adduct="M-H", Tolerance= 0.25, title="",Neighbour = 0,Creat_new_file=TRUE,color="black",interpolate =FALSE){
#x11()
#dev.copy(png,paste(dir,"/",png_filename,'.png',sep=""))
Tolerance<-round(Tolerance,digits = 5)
pngfillewrite<-paste(WKdir,"/",png_filename,'.png',sep="")
#pngfillewrite<-nice_file_create(pngfillewrite)
png(paste(WKdir,"/","temp.png",sep=""), bg = "transparent")
try(plotMsiSlice(imagefile ,mz , tolerance=Tolerance, legend=FALSE,colRamp=colorRamp(c("black", "blue", "green", "yellow", "red","#FF00FF","white")),interpolate =interpolate ),silent = TRUE)
dev.off()
try(pngfile<-image_read(paste(WKdir,"/","temp.png",sep="")))
res <- try(pngfile<-image_trim(pngfile),silent = TRUE)
if (class(res) != "try-error"){
a<-image_attributes(pngfile)
pixelx<-as.numeric(str_split(a[a$property=="png:IHDR.width,height","value"],", ")[[1]][1])
pixely<-as.numeric(str_split(a[a$property=="png:IHDR.width,height","value"],", ")[[1]][2])
try(resolutionx<-pixelx/imagefile[[1]]@metaData[["imaging"]][["size"]]["x"])
try(resolutiony<-pixely/imagefile[[1]]@metaData[["imaging"]][["size"]]["y"])
try(resolution<-sqrt(as.numeric((resolutionx + resolutiony) / 2)))
try(kern<-Kernel_convolution(Neighbour,resolution))
if (title=="D-mode"){
pngfile<-image_border(pngfile, "transparent", "0x30")
pngfile<-image_annotate(pngfile,paste(png_filename,"D-mode","\n",adduct,sprintf("%g",mz),"±",Tolerance),gravity = "south",size = 9)
pngfile<-image_trim(pngfile)
res2 <- try(pngfileoriginal<-image_read(pngfillewrite),silent = TRUE)
if ('&'(class(res2) != "try-error",Creat_new_file==FALSE)){
pngfile<-image_append(c(pngfileoriginal,pngfile))
image_write(pngfile,pngfillewrite)}
if ('&'(class(res2) != "try-error",Creat_new_file==TRUE)){
image_write(pngfile,pngfillewrite)}
if (class(res2) == "try-error"){
image_write(pngfile,pngfillewrite)}
}else{
pngfile<-image_convolve(pngfile, kern)
pngfile<-image_border(pngfile, "transparent", "20x20")
pngfile<-image_annotate(pngfile,paste(png_filename,"\n", adduct,sprintf("%g",mz),"±",Tolerance),gravity = "south",size = 10,color = color)
pngfile<-image_trim(pngfile)
image_write(pngfile,pngfillewrite)
}
}
file.remove(paste(WKdir,"/","temp.png",sep=""))
}
Plot_Ion_image_Png2<- function(WKdir, imagefile, png_filename, mz,adduct="M-H", Tolerance= 0.25, title="",Creat_new_file=TRUE){
#x11()
pngfillewrite<-paste(WKdir,"/",png_filename,'.png',sep="")
#pngfillewrite<-nice_file_create(pngfillewrite)
#dev.copy(png,paste(dir,"/",png_filename,'.png',sep=""))
png(paste(WKdir,"/","temp.png",sep=""))
try(plotMsiSlice(imagefile ,mz , tolerance=Tolerance,colRamp=colorRamp(c("black", "blue", "green", "yellow", "red")), legend=FALSE),silent = TRUE)
dev.off()
try(pngfile<-image_read(paste(WKdir,"/","temp.png",sep="")))
res <- try(pngfile<-image_trim(pngfile),silent = TRUE)
if (class(res) != "try-error"){
if (title=="D-mode"){
pngfile<-image_annotate(pngfile,paste(png_filename,"D-mode","\n",adduct,sprintf("%g",mz),"±",Tolerance),size = 9)
res2 <- try(pngfileoriginal<-image_read(pngfillewrite),silent = TRUE)
if ('&'(class(res2) != "try-error",Creat_new_file==FALSE)){
pngfile<-image_append(c(pngfileoriginal,pngfile))
image_write(pngfile,pngfillewrite)}
if ('&'(class(res2) != "try-error",Creat_new_file==TRUE)){
image_write(pngfile,pngfillewrite)}
if (class(res2) == "try-error"){
image_write(pngfile,pngfillewrite)}
}else{
pngfile<-image_annotate(pngfile,paste(png_filename,"\n", adduct,sprintf("%g",mz),"±",Tolerance),size = 9)
image_write(pngfile,pngfillewrite)
}
}
file.remove(paste(WKdir,"/","temp.png",sep=""))
}
#wrap the getMolecule function
##prepare SMPDB
#prepare_SMPDB<- function(DBpath="DB"){
#SMPDB <- list.files(path=paste(DBpath,"/metapath",sep=""),full.names = TRUE) %>% lapply(read.csv) %>% bind_rows
#write.csv(SMPDB, file=paste(DBpath,"/SMPdb.csv",sep=","))
#if (exists("SMPDB")==FALSE){SMPDB <- read.csv(file="DB/SMPdb.csv", as.is = TRUE)}
#SMPDBFREQ<-as.data.frame(table(SMPDB[,"SMPDB.ID"]))
#write.csv(SMPDBFREQ, file=(paste(DBpath,"/SMPdb_ID_FREQ.csv",sep="")))
#DBpath="DB"
#SMPDB <- read.csv(file=paste(DBpath,"/SMPdb.csv",sep=""), as.is = TRUE)
#SMPDBFREQ <- read.csv(file=paste(DBpath,"/SMPdb_ID_FREQ.csv",sep="") ,col.names= c("","SMPDB.ID","Freq"),as.is = TRUE)
#pathwayDB <- read.csv(file=paste(DBpath,"/smpdb_pathways.csv",sep=""), as.is = TRUE)
#R Python interface
#install.packages("reticulate")
#install.packages("R.matlab")
#}
quiet <- function(x) {
sink(tempfile())
on.exit(sink())
invisible(force(x))
}
Cardinal_old<- function(){
workdir <- WorkingDir()
datafile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.imzML$")
listfile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.csv$")
datafile<-gsub(".imzML", "", datafile)
name <-gsub(paste(workdir,"\\/",sep=""),"",datafile[3])
#imdata <- readImzML(name, folder, as="MSImagingExperiment")
imdata <- readImzML(name, folder, attach.only=TRUE, as="MSImagingExperiment")
summarize(data, sum, .by="pixel")
tmp <- imdata %>%
smoothSignal() %>%
reduceBaseline() %>%
peakPick() %>%
peakFilter() %>%
select(x == 1, y == 1)
process(plot=TRUE,
par=list(layout=c(1,3)),
BPPARAM=SerialParam())
plot(data, pixel=1)
plot(data, coord=list(x=2, y=2))
pattern <- factor(c(0, 0, 2, 2, 0, 0, 0, 0, 0, 0, 2, 2, 0,
0, 0, 0, 0, 0, 0, 1, 2, 2, 0, 0, 0, 0, 0, 2, 1, 1, 2,
2, 0, 0, 0, 0, 0, 1, 2, 2, 2, 2, 0, 0, 0, 0, 1, 2, 2,
2, 2, 2, 0, 0, 0, 0, 2, 2, 2, 2, 2, 2, 2, 0, 0, 0, 2,
2, 0, 0, 0, 0, 0, 0, 2, 2, 0, 0, 0, 0, 0),
levels=c(0,1,2), labels=c("blue", "black", "red"))
tmp <- imdata %>%
smoothSignal() %>%
reduceBaseline() %>%
peakPick() %>%
peakFilter() %>%
select(x == 1, y == 1) %>%
process(plot=TRUE,
par=list(layout=c(1,3)),
BPPARAM=SerialParam())
}
intensity.colors_customize <- function(n = 100, alpha = 1) {
col2 <- rainbow(3*n, alpha=alpha)[(2*n):1]
f <- colorRamp(colors=c("black", "blue", "green", "yellow", "red","#FF00FF","white"))
alpha <- col2rgb(col2, alpha=TRUE)[[4]]
col1 <- sapply(seq(from=0, to=1, length.out=n), function(i) do.call(rgb,
c(as.list(f(i)), maxColorValue=255, alpha=alpha)))
col1
}
Pathway_overview_graphite<-function(){
p_load(graphite,graph )
p_load(org.Hs.eg.db)
p_load(Rgraphviz)
humanReactome <- graphite::pathways("hsapiens", "reactome")
humanmeatbolome <- graphite::pathways("hsapiens", "smpdb")
metab_url <-
url("https://romualdi.bio.unipd.it/wp-uploads/2018/04/Terunuma_metabolite_expr.txt")
mexpr <- read.table(metab_url, header = TRUE, sep = "\t", row.names = NULL,
stringsAsFactors = FALSE)
p <- humanReactome[["Glycolysis"]]
head(nodes(p))
head(edges(p))
head(nodes(p), which = "mixed")
head(edges(p), which = "mixed")
pathwayDatabases()
g <- pathwayGraph(p)
g <- pathwayGraph(p, which = "mixed")
pSymbol <- convertIdentifiers(p, "SYMBOL")
head(nodes(pSymbol))
reactomeSymbol <- convertIdentifiers(humanReactome[1:5], "SYMBOL")
cytoscapePlot(convertIdentifiers(reactome$`Unwinding of DNA`, "symbol"), which = "mixed")
meta_g <- pathwayGraph(p, which = "metabolites")
plot(g)
}
Cardinal_cluster_classification<-function(){
set.seed(1)
skm <- spatialKMeans(imdata, r=2, k=4, method="adaptive")
plot(skm, col=c("pink", "blue", "red","orange","navyblue"), type=c('p','h'), key=FALSE)
image(skm, col=c("pink", "blue", "red","orange","navyblue"), key=FALSE)
set.seed(1)
ssc <- spatialShrunkenCentroids(imdata, r=2, k=5, s=3, method="gaussian")
plot(ssc, col=c("pink", "blue", "red","orange","navyblue"), type=c('p','h'), key=FALSE)
image(ssc, col=c("pink", "blue", "red","orange","navyblue"), key=FALSE)
}
MASHUOA/HiTMaP documentation built on May 1, 2024, 8:26 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Defines functions Cardinal_cluster_classification Pathway_overview_graphite intensity.colors_customize Cardinal_old quiet Plot_Ion_image_Png2 Plot_Ion_image_Png3 Plot_Ion_image_Png Kernel_convolution nice_file_create fetch_ID_from_DB Whole_Picture_Cluster_mode META_pathway_Whole_Picture_discovery_mode META_pathway_Whole_Picture META_pathway_ion_image_cardinal META_pathway_ion_image META_pathway_search WorkingDir simple_ion_image Generate_HMDBIDS cluster_image_cardinal intensity.colors_customize_mono metwork_cardinal metwork_maldiquant Cardinal_utilities Generatespectrum_workflow_maldiquant Cardinal run_spatial_quant_workflow imaging_default_parameter RUN_imagine_workflow
# code for testing purpose
RUN_imagine_workflow<-function(){
#Spatial_Quant(ppm=15,adducts=c("M-H","M+Cl"),Quant_list="lipid candidates.csv",cal.mz=T)
#Spatial_Quant(ppm=15,adducts=c("M-H","M+Cl"),Quant_list="lipid candidates.csv",cal.mz=T,Spectrum_feature_summary = T,Region_feature_summary = T,PMF_analysis = F)
imaging_identification(
#==============Choose the imzml raw data file(s) to process
datafile=tk_choose.files(filter = matrix(c( "imzml file", ".imzML",
"Text", ".txt", "All files", "*"),
3, 2, byrow = TRUE),caption = "Choose single or multiple file(s) for analysis"),
threshold=0.05,
ppm=20,
Digestion_site="[G]",
missedCleavages=0:2,
Fastadatabase="murine_matrisome.fasta",
adducts=c("M+H","M+NH4","M+Na"),
#==============TRUE if you want to perform PMF analysis
PMF_analysis=T,
#==============TRUE if you want to generate peptide summary in the Summary folder
Protein_feature_summary=T,
Peptide_feature_summary=T,
#==============TRUE if you want to plot peptide in the Ion images folder, make sure there's imzml file in the folder
plot_ion_image=F,
#==============Set a number if you want a parallel processing
parallel=8,
#==============Set a number (1 to maximum pixels in the data file) if you want to dig more peaks in the raw data
spectra_segments_per_file=5,
#==============True if you want the mean spectrum generated based on spatial spectrum profile
spatialKMeans=T,
Smooth_range=1
)
}
imaging_default_parameter<-function(){
datafile=tk_choose.files(filter = matrix(c( "imzml file", ".imzML",
"Text", ".txt", "All files", "*"),
3, 2, byrow = TRUE),caption = "Choose single or multiple file(s) for analysis")
threshold=0.05
ppm=10
Digestion_site="[G]"
missedCleavages=0:2
Fastadatabase="murine_matrisome.fasta"
adducts=c("M+H","M+NH4","M+Na")
#==============TRUE if you want to perform PMF analysis
PMF_analysis=T
#==============TRUE if you want to generate peptide summary in the Summary folder
Protein_feature_summary=T
Peptide_feature_summary=T
#==============TRUE if you want to plot peptide in the Ion images folder, make sure there's imzml file in the folder
plot_ion_image=T
#==============Set a number if you want a parallel processing
parallel=8
#==============Set a number (1 to maximum pixels in the data file) if you want to dig more peaks in the raw data
spectra_segments_per_file=5
spatialKMeans=T
plot_cluster_image=T
}
run_spatial_quant_workflow<-function(){
imagine_Spatial_Quant(datafile=tk_choose.files(filter = matrix(c( "imzml file", ".imzML",
"Text", ".txt", "All files", "*"),
3, 2, byrow = TRUE),caption = "Choose single or multiple file(s) for analysis"),
threshold=0.00,
ppm=5,
#Fastadatabase="murine_matrisome.fasta",
Quant_list="Metabolites of Interest.csv",
adducts=c("M-H","M+Cl"),
cal.mz=F,
mzlist_bypass=T,
#==============TRUE if you want to plot protein PMF result
PMF_analysis=TRUE,
#==============TRUE if you want to generate peptide summary in the Summary folder
Protein_feature_summary=F,
Peptide_feature_summary=F,
#PLOT_PMF_Protein=FALSE,
#==============TRUE if you want to plot peptide in the Ion images folder, make sure there's imzml file in the folder
plot_ion_image=FALSE,
#==============Set a number if you want a parallel processing
parallel=detectCores()/2,
#==============Set a number (1 to maximum pixels in the data file) if you want to dig more peaks in the raw data
spectra_segments_per_file=4,
spatialKMeans=F,
Smooth_range=1,
Virtual_segmentation=T,
Virtual_segmentation_rankfile="Z:\\George skyline results\\maldiimaging\\Maldi_imaging - Copy\\radius_rank_bovin.csv",
Spectrum_feature_summary=F,
Region_feature_summary=F,
Region_feature_analysis=F,
Cluster_level=="Low")
#workflow for aged human lens
imagine_Spatial_Quant(
#==============Choose the imzml raw data file(s) to process make sure the fasta file in the same folder
datafile=tk_choose.files(filter = matrix(c( "imzml file", ".imzML",
"Text", ".txt", "All files", "*"),
3, 2, byrow = TRUE),caption = "Choose single or multiple file(s) for analysis"),
threshold=0.00,
ppm=20,
#Fastadatabase="murine_matrisome.fasta",
Quant_list="lipid candidates.csv",
adducts=c("M-H","M+Cl"),
cal.mz=T,
mzlist_bypass=F,
#==============TRUE if you want to plot protein PMF result
PMF_analysis=TRUE,
#==============TRUE if you want to generate protein summary in the Summary folder
Protein_feature_summary=T,
#==============TRUE if you want to generate protein cluster image in the Summary folder
plot_cluster_image=T,
Peptide_feature_summary=T,
#PLOT_PMF_Protein=FALSE,
#==============TRUE if you want to plot peptide in the Ion images folder, make sure there's imzml file in the folder
plot_ion_image=FALSE,
#==============Set a number if you want a parallel processing
parallel=detectCores(),
#==============Set a number (1 to maximum pixels in the data file) if you want to dig more peaks in the raw data
spectra_segments_per_file=5,
spatialKMeans=F,
Smooth_range=1,
Virtual_segmentation=T,
Virtual_segmentation_rankfile="Z:\\George skyline results\\maldiimaging\\Maldi_imaging - Copy\\radius_rank.csv",
Spectrum_feature_summary=T,
Region_feature_summary=T,
Region_feature_analysis=T,
plot_each_metabolites=F,
Cluster_level="High",
Region_feature_analysis_bar_plot=F,
norm_datafiles=T,
norm_Type="Median"
)
imagine_Spatial_Quant(
#==============Choose the imzml raw data file(s) to process make sure the fasta file in the same folder
datafile=tk_choose.files(filter = matrix(c( "imzml file", ".imzML",
"Text", ".txt", "All files", "*"),
3, 2, byrow = TRUE),caption = "Choose single or multiple file(s) for analysis"),
threshold=0.00,
ppm=25,
#Fastadatabase="murine_matrisome.fasta",
Quant_list="lipid blast candidates.csv",
adducts=c("M-H","M+Cl"),
cal.mz=F,
mzlist_bypass=F,
#==============TRUE if you want to plot protein PMF result
PMF_analysis=TRUE,
#==============TRUE if you want to generate protein summary in the Summary folder
Protein_feature_summary=T,
#==============TRUE if you want to generate protein cluster image in the Summary folder
plot_cluster_image=T,
Peptide_feature_summary=T,
#PLOT_PMF_Protein=FALSE,
#==============TRUE if you want to plot peptide in the Ion images folder, make sure there's imzml file in the folder
plot_ion_image=FALSE,
#==============Set a number if you want a parallel processing
parallel=detectCores(),
#==============Set a number (1 to maximum pixels in the data file) if you want to dig more peaks in the raw data
spectra_segments_per_file=5,
spatialKMeans=F,
Smooth_range=1,
Virtual_segmentation=T,
Virtual_segmentation_rankfile="Z:\\George skyline results\\maldiimaging\\Maldi_imaging - Copy\\radius_rank.csv",
Spectrum_feature_summary=T,
Region_feature_summary=T,
Region_feature_analysis=T,
plot_each_metabolites=F,
Cluster_level="High",
Region_feature_analysis_bar_plot=F,
norm_datafiles=T
)
}
Cardinal<- function(){
workdir <- WorkingDir()
datafile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.imzML$")
listfile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.csv$")
datafile<-gsub(".imzML", "", datafile)
name <-gsub(paste(workdir,"\\/",sep=""),"",datafile[3])
#imdata <- readImzML(name, folder, as="MSImagingExperiment")
imdata <- readImzML(name, folder, attach.only=TRUE, as="MSImagingExperiment")
summarize(data, sum, .by="pixel")
tmp <- imdata %>%
smoothSignal() %>%
reduceBaseline() %>%
peakPick() %>%
peakFilter() %>%
select(x == 1, y == 1)
process(plot=TRUE,
par=list(layout=c(1,3)),
BPPARAM=SerialParam())
plot(data, pixel=1)
plot(data, coord=list(x=2, y=2))
pattern <- factor(c(0, 0, 2, 2, 0, 0, 0, 0, 0, 0, 2, 2, 0,
0, 0, 0, 0, 0, 0, 1, 2, 2, 0, 0, 0, 0, 0, 2, 1, 1, 2,
2, 0, 0, 0, 0, 0, 1, 2, 2, 2, 2, 0, 0, 0, 0, 1, 2, 2,
2, 2, 2, 0, 0, 0, 0, 2, 2, 2, 2, 2, 2, 2, 0, 0, 0, 2,
2, 0, 0, 0, 0, 0, 0, 2, 2, 0, 0, 0, 0, 0),
levels=c(0,1,2), labels=c("blue", "black", "red"))
tmp <- imdata %>%
smoothSignal() %>%
reduceBaseline() %>%
peakPick() %>%
peakFilter() %>%
select(x == 1, y == 1) %>%
process(plot=TRUE,
par=list(layout=c(1,3)),
BPPARAM=SerialParam())
}
Generatespectrum_workflow_maldiquant<-function(){
table(sapply(MALDI_IMAGE, length))
plot(MALDI_IMAGE[[1]])
workdir <- WorkingDir()
datafile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.imzML$")
datafile<-gsub(".imzML", "", datafile)
MALDI_IMAGE <- quiet(importImzMl(paste(file.path( datafile[i]),".imzml",sep="")))
#use the square root transformation to simplify graphical visualization and to overcome the potential dependency of the variance from the mean.
spectra <- transformIntensity(MALDI_IMAGE, method="sqrt")
#use a 21 point Savitzky-Golay-Filter (Savitzky and Golay, 1964) to smooth the spectra.
spectra <- smoothIntensity(spectra, method="SavitzkyGolay", halfWindowSize=10)
#Before we correct the baseline we visualize it. Here we use the SNIP algorithm (Ryan et al., 1988).
baseline <- estimateBaseline(spectra[[16]], method="SNIP", iterations=100)
plot(spectra[[16]])
lines(baseline, col="red", lwd=2)
#If we are satisfied with our estimated baseline we remove it.
spectra <- removeBaseline(spectra, method="SNIP",iterations=100)
plot(spectra[[1]])
#spectra <- calibrateIntensity(spectra, method="TIC")
#do the alignment if there's a huge mass shift between each spectrum
spectra <- alignSpectra(spectra,halfWindowSize=50,SNR=3,tolerance=0.05,warpingMethod="lowess")
samples <- factor(sapply(spectra, function(x)metaData(x)$sampleName))
avgSpectra <- averageMassSpectra(spectra, labels=samples, method="mean")
noise <- estimateNoise(avgSpectra[[1]])
plot(avgSpectra[[1]], xlim=c(500, 4000), ylim=c(0, 0.002))
lines(noise, col="red")
lines(noise[,1], noise[, 2]*2, col="blue")
peaks <- detectPeaks(avgSpectra, method="MAD",halfWindowSize=20, SNR=2)
plot(avgSpectra[[1]], xlim=c(4000, 5000), ylim=c(0, 0.002))
points(peaks[[1]], col="red", pch=4)
peaks <- binPeaks(peaks, tolerance=0.002)
peaks <- filterPeaks(peaks, minFrequency=0.25)
}
Cardinal_utilities<- function(){
folder <- WorkingDir()
datafile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.imzML$")
listfile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.csv$")
datafile<-gsub(".imzML", "", datafile)
name <-gsub(paste(folder,"\\/",sep=""),"",datafile[1])
#imdata <- readImzML(name, folder, as="MSImagingExperiment")
imdata <- readImzML(name, folder, attach.only=TRUE, as="MSImagingExperiment")
summarize(data, sum, .by="pixel")
tmp <- data %>%
smoothSignal() %>%
reduceBaseline() %>%
peakPick() %>%
peakFilter() %>%
select(x == 1, y == 1)
process(plot=TRUE,
par=list(layout=c(1,3)),
BPPARAM=SerialParam())
plot(data, pixel=1)
plot(data, coord=list(x=2, y=2))
pattern <- factor(c(0, 0, 2, 2, 0, 0, 0, 0, 0, 0, 2, 2, 0,
0, 0, 0, 0, 0, 0, 1, 2, 2, 0, 0, 0, 0, 0, 2, 1, 1, 2,
2, 0, 0, 0, 0, 0, 1, 2, 2, 2, 2, 0, 0, 0, 0, 1, 2, 2,
2, 2, 2, 0, 0, 0, 0, 2, 2, 2, 2, 2, 2, 2, 0, 0, 0, 2,
2, 0, 0, 0, 0, 0, 0, 2, 2, 0, 0, 0, 0, 0),
levels=c(0,1,2), labels=c("blue", "black", "red"))
tmp <- imdata %>%
smoothSignal() %>%
reduceBaseline() %>%
peakPick() %>%
peakFilter() %>%
select(x == 1, y == 1) %>%
process(plot=TRUE,
par=list(layout=c(1,3)),
BPPARAM=SerialParam())
}
####load pathway database
#Map_SMP_MP <- read.csv("DB//smpdb_pathways.csv",as.is = TRUE)
####load raw data annotation
####load adducts data
#adducts<-read.csv("DB//adducts.csv",as.is = TRUE)
####load SMPDB data
#SMPDB <- read.csv(file="DB/SMPdb.csv",stringAsFactors=FALSE,as.is = TRUE)
#SMPDBFREQ <- read.csv(file="DB/SMPdb_ID_FREQ.csv",col.names= c("","SMPDB.ID","Freq"),as.is = TRUE)
##work dir
##load raw data file
#MALDI_IMAGE <- importImzMl(file.path( "NEDC.imzML"))
######Example steps
metwork_maldiquant<- function(workdir= tk_choose.dir(caption = "Select working directory")){
datafile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.imzML$")
datafile<-gsub(".imzML", "", datafile)
setwd("Z:\\George skyline results\\maldiimaging\\Maldi_imaging - Copy")
HMDBids <- read.csv("combine.csv", fill = TRUE,skip=2,as.is = TRUE)
SMPDB <- read.csv(file="DB/SMPdb.csv",as.is = TRUE)
for (i in 1: length(datafile)){
if (dir.exists(datafile[i])==FALSE){dir.create(datafile[i])}
HMDBDATA <- Generate_HMDBIDS(HMDBids)
SMPLIST <- META_pathway_search(HMDBDATA, SMPDB, datafile[i], PWtype="Metabolic")
SMPRESULT <-read.csv(file=paste(datafile[i],sep="","/","SMPRESULT.csv"),as.is = TRUE)
options(warn=-1)
MALDI_IMAGE <- quiet(importImzMl(paste(file.path( datafile[i]),".imzml",sep="")))
options(warn=0)
META_pathway_ion_image(SMPLIST, SMPRESULT, datafile[i] ,MALDI_IMAGE, Image_Type="",PWpScore=0.26,plotTolerance=0.125 )
META_pathway_Whole_Picture(SMPLIST, SMPRESULT, datafile[i], PWpScore=0.26,pathway_mode="")
META_pathway_Whole_Picture_discovery_mode(SMPLIST, SMPRESULT, datafile[i], SMPDB, adductslist=c("M-H","M+Cl"), PWpScore=0.26)
}
}
metwork_cardinal<- function(workdir= tk_choose.dir(caption = "Select working directory")){
datafile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.imzML$")
datafile<-gsub(".imzML", "", datafile)
setwd("Z:\\George skyline results\\maldiimaging\\Maldi_imaging - Copy")
HMDBids <- read.csv("combine.csv", fill = TRUE,skip=2,as.is = TRUE)
SMPDB <- read.csv(file="DB/SMPdb.csv",as.is = TRUE)
for (i in 1: length(datafile)){
if (dir.exists(datafile[i])==FALSE){dir.create(datafile[i])}
HMDBDATA <- Generate_HMDBIDS(HMDBids)
SMPLIST <- META_pathway_search(HMDBDATA, SMPDB, datafile[i], PWtype="Metabolic")
SMPRESULT <-read.csv(file=paste(datafile[i],sep="","/","SMPRESULT.csv"),as.is = TRUE)
options(warn=-1)
imdata <- quiet(Load_Cardinal_imaging(paste(file.path( datafile[i]),sep=""),resolution = 20))
options(warn=0)
META_pathway_ion_image_cardinal(SMPLIST, SMPRESULT, datafile[i] ,imdata, Image_Type="",PWpScore=0.26,plotTolerance=20 ,)
filename=gsub(dirname(datafile[i]),"",datafile[i])
for (SMPDB.ID in SMPRESULT$SMPDB.ID[SMPRESULT$pScore>=PWpScore]){
cluster_image_cardinal(SMPDB.ID,imdata,SMPLIST,SMPRESULT$Name[SMPRESULT$SMPDB.ID==SMPDB.ID],ppm=ppm)
#cluster_image_cardinal(SMPDB.ID,imdata_istd,SMPLIST,SMPRESULT$Name[SMPRESULT$SMPDB.ID==SMPDB.ID])
}
#lapply(SMPRESULT$SMPDB.ID[SMPRESULT$pScore>=PWpScore],cluster_image_cardinal,imdata,SMPLIST)
META_pathway_Whole_Picture(SMPLIST, SMPRESULT, datafile[i], PWpScore=0.26,pathway_mode="")
META_pathway_Whole_Picture_discovery_mode(SMPLIST, SMPRESULT, datafile[i], SMPDB, adductslist=c("M-H","M+Cl"), PWpScore=0.26)
}
}
##get HMDB IDs list for pathway search
intensity.colors_customize_mono <- function(n = 100, alpha = 1,colorstart="transparent",colorend="red") {
col1 <- colorRamp(colors=c(colorstart, colorend),interpolate = "linear")
col1
}
cluster_image_cardinal<-function(clusterID,imdata,SMPLIST,clustername,ppm=20){
library(RColorBrewer)
outputpng=paste(getwd(),"/",clustername,"_cluster_plot",'.png',sep="")
candidate=SMPLIST[SMPLIST$SMPDB.ID==clusterID,]
candidateunique=unique(candidate[,"mz"])
mycol <- RColorBrewer::brewer.pal(length(candidateunique),"Set1")
png(outputpng,width = 720,height = 720, bg = "transparent")
image(imdata, mz=candidateunique,
col=mycol,
#contrast.enhance = contrast.enhance,
smooth.image = smooth.image ,
superpose=TRUE,normalize.image="linear")
dev.off()
pngfile<-image_read(outputpng)
pngfile<-image_border(pngfile, "transparent", "30x30")
pngfile<-image_annotate(pngfile,paste(clustername),gravity = "south",size = 40,color = "white")
pngfile<-image_trim(pngfile)
image_write(pngfile,outputpng)
for (i in 1:length(candidateunique)){
outputpngsum=paste(getwd(),"/",clustername,"_cluster_plot_",candidateunique[i],'.png',sep="")
png(outputpngsum,width = 720,height = 720, bg = "transparent")
par(oma=c(0, 0, 0, 0),tcl = NA,mar=c(0, 0, 1, 1),mfrow = c(1, 1),
bty="n",pty="s",xaxt="n",
yaxt="n",
no.readonly = TRUE,ann=FALSE)
#image(imdata, mz=candidateunique[i], col=mycol[i], superpose=F,normalize.image="linear")
col.regions <- gradient.colors(100, start="transparent", end=mycol[i])
image(imdata, mz=candidateunique[i],
contrast.enhance=contrast.enhance,
smooth.image = smooth.image ,
col.regions=col.regions,normalize.image="linear")
dev.off()
}
}
#plot cluster image all in one file
Generate_HMDBIDS <- function(HMDBids){
HMDBDATA<-NULL
HMDBSearchlist<-str_split(paste(HMDBids$moleculeIds,collapse=", "),pattern=", ")[[1]]
Formular<-str_split(paste(HMDBids$moleculeNames,collapse=", "),pattern=", ")[[1]]
HMDBDATA$HMDB.ID<-HMDBSearchlist
HMDBDATA$moleculeNames<-Formular
MZLIST <- NULL
adduct <- NULL
options(warn=-1)
for (i in 1: length(HMDBSearchlist)){
MZLIST[i]<-as.numeric(HMDBids[grep(HMDBSearchlist[i],HMDBids[,"moleculeIds"]),"mz"])
adduct[i]<-as.character(HMDBids[grep(HMDBSearchlist[i],HMDBids[,"moleculeIds"]),"adduct"])
}
HMDBDATA$mz<-MZLIST
HMDBDATA$adduct<-adduct
options(warn=0)
HMDBDATA<-as.data.frame(HMDBDATA)
HMDBDATA<-unique(HMDBDATA)
return(HMDBDATA)
}
##simple_ion_image: define the work dir read all csv and enerate the ion image respectively
simple_ion_image<-function(workdir=WorkingDir(),
Image_Type="",
plotTolerance=5 ,
creat_new_file=TRUE,
Denoising=0,
color = "black",
interpolate =FALSE,
...){
datafile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.imzML$")
listfile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.csv$")
masslist<-NULL
for (z in 1:length(datafile)){
datafile[z]<-gsub(".imzML", "", datafile[z])
options(warn=-1)
MALDI_IMAGE <- importImzMl(paste(file.path(datafile[z]),".imzML",sep=""))
options(warn=0)
if (dir.exists(datafile[z])==FALSE){dir.create(datafile[z])}
for (i in 1:length(listfile)){
foldername<-str_remove(listfile[i],paste(workdir,"/",sep=""))
foldername<-str_remove(foldername,".csv")
outputDir<-paste(datafile[z],"/",foldername,sep="")
if (dir.exists(outputDir)==FALSE){dir.create(outputDir)}
masslist<-read.csv(listfile[i],header = TRUE)
for (j in 1:length(masslist$moleculeNames)){
if (is.na(masslist$mz[j])==FALSE){
for (k in 1:length(str_split(masslist$moleculeNames[j],", ")[[1]])){
Plot_Ion_image_Png3(outputDir, MALDI_IMAGE, windows_filename(str_split(masslist$moleculeNames[j],", ")[[1]][k]), masslist$mz[j],masslist$adduct[j], Tolerance= plotTolerance*masslist$mz[j]/1000000, title="",Creat_new_file=TRUE,Neighbour = Denoising,color = color,interpolate =interpolate)
}}}
}
}
}
WorkingDir <-function(Mainfolder = tk_choose.dir(caption = "Select working directory")){
Mainfolder <- Mainfolder
setwd(Mainfolder)
Mainfolder
}
##Pathway search return search result
META_pathway_search <- function(HMDBDATA, SMPDB, WKdir, PWtype="Metabolic", ...) {
SMPLIST<- NULL
SMPMZlist<-NULL
for (i in 1: length(HMDBDATA[,"HMDB.ID"])){
SMPquery<-NULL
SMPquery<-SMPDB[`&`(SMPDB[,"Pathway.Subject"]==PWtype,SMPDB[,"HMDB.ID"] == HMDBDATA[,"HMDB.ID"][i]),"SMPDB.ID"]
SMPLIST$SMPDB.ID<-c(SMPLIST$SMPDB.ID,SMPquery)
SMPLIST$HMDB.ID<-c(SMPLIST$HMDB.ID, rep(as.character(HMDBDATA[,"HMDB.ID"][i]),length(SMPquery)))
SMPLIST$mz<-c(SMPLIST$mz, rep(HMDBDATA[,"mz"][i],length(SMPquery)))
SMPLIST$moleculeNames<-c(SMPLIST$moleculeNames, rep(as.character(HMDBDATA[,"moleculeNames"][i]),length(SMPquery)))
SMPLIST$KEGG.ID<-c(SMPLIST$KEGG.ID, rep(SMPDB[SMPDB[,"HMDB.ID"]==HMDBDATA[,"HMDB.ID"][i],"KEGG.ID"][1],length(SMPquery)))
SMPLIST$Metabolite.ID<-c(SMPLIST$Metabolite.ID, rep(SMPDB[SMPDB[,"HMDB.ID"]==HMDBDATA[,"HMDB.ID"][i],"Metabolite.ID"][1],length(SMPquery)))
SMPLIST$adduct<-c(SMPLIST$adduct, rep(as.character(HMDBDATA[HMDBDATA[,"HMDB.ID"]==HMDBDATA[,"HMDB.ID"][i],"adduct"][1]),length(SMPquery)))
}
if (dir.exists(paste(WKdir,sep="","/"))==FALSE){dir.create(paste(WKdir,sep="","/"))}
SMPLIST<-as.data.frame(SMPLIST)
write.csv(SMPLIST,paste(WKdir,sep="","/","SMPqueryRESULT.csv"),row.names=FALSE)
SMPRESULT<-as.data.frame(table(SMPLIST$SMPDB.ID))
write.csv(SMPRESULT,paste(WKdir,sep="","/","SMPRESULT.csv"),row.names=FALSE)
SMPRESULT<-read.csv(file=paste(WKdir,sep="","/","SMPRESULT.csv"),col.names= c("SMPDB.ID","Freq"))
SMPRESULTfreq<-NULL
for (i in 1: length(SMPRESULT$SMPDB.ID)){
SMPRESULTfreq<-c(SMPRESULTfreq,SMPDBFREQ[SMPDBFREQ[,"SMPDB.ID"]==as.character(SMPRESULT[i,"SMPDB.ID"]),"Freq"])
SMPRESULT$Name[i]<-as.character(Map_SMP_MP[Map_SMP_MP[,"SMPDB.ID"] %in% SMPRESULT[,"SMPDB.ID"][i],"Name"])
}
SMPRESULT$DBfreq<-SMPRESULTfreq
SMPRESULT$pScore<-SMPRESULT$Freq/SMPRESULT$DBfreq
write.csv(SMPRESULT,paste(WKdir,sep="","/","SMPRESULT.csv"),row.names=FALSE)
return(SMPLIST)
}
###Ion Image manipulation
META_pathway_ion_image<- function(SMPLIST, SMPRESULT, datafile,MALDI_IMAGE, Image_Type="",PWpScore=0.26,plotTolerance=0.125 ,creat_new_file=TRUE,...) {
if (dir.exists(datafile)==FALSE){dir.create(datafile)}
for (j in 1:length(SMPRESULT$SMPDB.ID)){
if (SMPRESULT$pScore[j]>=PWpScore){
WKdir<-paste(datafile,"/",SMPRESULT$SMPDB.ID[j],sep="")
PWSUBLIST<-SMPLIST[SMPLIST$SMPDB.ID==SMPRESULT$SMPDB.ID[j],]
if (dir.exists(WKdir)==FALSE){dir.create(WKdir)}
for (i in 1:length(PWSUBLIST$mz)){
Plot_Ion_image_Png(WKdir,MALDI_IMAGE,PWSUBLIST$moleculeNames[i],PWSUBLIST$mz[i],PWSUBLIST$adduct[i],Tolerance= plotTolerance, title=Image_Type,Creat_new_file=creat_new_file)
}}
}}
META_pathway_ion_image_cardinal<- function(SMPLIST, SMPRESULT, datafile,MALDI_IMAGE, Image_Type="",PWpScore=0.26,plotTolerance=0.125 ,creat_new_file=TRUE,...) {
if (dir.exists(datafile)==FALSE){dir.create(datafile)}
for (j in 1:length(SMPRESULT$SMPDB.ID)){
if (SMPRESULT$pScore[j]>=PWpScore){
WKdir<-paste(datafile,"/",SMPRESULT$SMPDB.ID[j],sep="")
PWSUBLIST<-SMPLIST[SMPLIST$SMPDB.ID==SMPRESULT$SMPDB.ID[j],]
if (dir.exists(WKdir)==FALSE){dir.create(WKdir)}
for (i in 1:length(PWSUBLIST$mz)){
Plot_Ion_image_Png_Cardinal(WKdir,MALDI_IMAGE,PWSUBLIST$moleculeNames[i],PWSUBLIST$mz[i],PWSUBLIST$adduct[i],Tolerance= plotTolerance, title=Image_Type,Creat_new_file=creat_new_file,contrast.enhance = "suppression")
}}
}}
###DRAW META_pathway_Whole_Picture
META_pathway_Whole_Picture<- function(SMPLIST, SMPRESULT, datafile, PWpScore=0.26,pathway_mode="", ...) {
MPID <-NULL
if (dir.exists(datafile)==FALSE){dir.create(datafile)}
for (j in 1:length(SMPRESULT$SMPDB.ID)){
if (SMPRESULT$pScore[j]>=PWpScore){
WKdir<-paste(datafile,"/",SMPRESULT$SMPDB.ID[j],sep="")
if (exists("Map_SMP_MP")==FALSE){Map_SMP_MP <- read.csv("DB/smpdb_pathways.csv")}
MPID <- sprintf("PW%06d", match(SMPRESULT$SMPDB.ID[j],Map_SMP_MP[,1]))
#PWsvg <- htmlParse(paste("DB/smpdb_svg/",MPID ,".svg",sep=""))
PWSUBLIST<-SMPLIST[SMPLIST$SMPDB.ID==SMPRESULT$SMPDB.ID[j],c("HMDB.ID","mz","moleculeNames","Metabolite.ID")]
PW_svg_line<- readLines(paste("DB/smpdb_svg/",MPID ,".svg",sep=""))
for (i in 1:length(PWSUBLIST$mz)){
pngfile<- NULL
res <- try(pngfile<-image_read(paste(WKdir,"/",PWSUBLIST$moleculeNames[i],'.png',sep="")),silent = TRUE)
if (class(res) != "try-error"){
Ion_image_svg_base64 <- base64Encode(readBin(paste(WKdir,"/",PWSUBLIST$moleculeNames[i],'.png',sep=""), "raw", file.info(paste(WKdir,"/",PWSUBLIST$moleculeNames[i],'.png',sep=""))[1, "size"]), "txt")
res<-try(Start_moldb_cache<-grep(paste('id=','"',PWSUBLIST$Metabolite.ID[i],'_moldb_cache',sep=""),as.character(PW_svg_line),value=FALSE),silent = TRUE)
res2<-try(End_moldb_cache<-grep("pointer-events=",as.character(PW_svg_line[Start_moldb_cache:length(PW_svg_line)]),value=FALSE)[1],silent = TRUE)
if ("&"(class(res) != "try-error",class(res2) != "try-error")){
PW_svg_line<-c(PW_svg_line[1:Start_moldb_cache], paste("<image xlink:href=\"data:image/png;base64,",Ion_image_svg_base64,aep=""),PW_svg_line[(Start_moldb_cache+End_moldb_cache-1):length(PW_svg_line)])
if (pathway_mode=="D-mode"){write(PW_svg_line,file=paste(datafile,"/",SMPRESULT$SMPDB.ID[j]," D-mode.svg",sep=""))}else
{write(PW_svg_line,file=paste(datafile,"/",SMPRESULT$SMPDB.ID[j],".svg",sep=""))}}
}
}
}
}
}
#Metapathway Image discovery_mode_list
META_pathway_Whole_Picture_discovery_mode<- function(SMPLIST, SMPRESULT, datafile, SMPDB, adductslist="M-H", PWpScore=0.26, ...) {
MPID <-NULL
SMPLIST_D_MODE<-SMPDB[`&`(SMPDB[,"SMPDB.ID"] %in% SMPRESULT[SMPRESULT[,"pScore"]>0.26,"SMPDB.ID"],!(SMPDB[,"HMDB.ID"] %in% unique(SMPLIST[,"HMDB.ID"]))),]
SMPLIST_D_MODE$mass<- t(as.data.frame(lapply(SMPLIST_D_MODE$Formula, getMonomass)))
SMPLIST_D_MODE$moleculeNames<-SMPLIST_D_MODE$Metabolite.Name
for (adduction in 1:length(adductslist)){
SMPLIST_D_MODE$adduct<-adductslist[adduction]
SMPLIST_D_MODE$mz <- SMPLIST_D_MODE$mass+ adducts[adducts[,"Name"]==adductslist[adduction],"Mass"]
if (adduction==1){
META_pathway_ion_image(SMPLIST_D_MODE, SMPRESULT, datafile,MALDI_IMAGE, Image_Type="D-mode", PWpScore=0.26,plotTolerance=0.125,creat_new_file=TRUE)
}else{
META_pathway_ion_image(SMPLIST_D_MODE, SMPRESULT, datafile,MALDI_IMAGE, Image_Type="D-mode", PWpScore=0.26,plotTolerance=0.125,creat_new_file=FALSE)
}
}
options(warn=-1)
SMPLIST_D_MODE_combine<-bind_rows(SMPLIST_D_MODE[,colnames(SMPLIST)],SMPLIST)
options(warn=0)
META_pathway_Whole_Picture(SMPLIST_D_MODE_combine, SMPRESULT, datafile, PWpScore=0.26,pathway_mode="D-mode")
}
#Cluster_mode_list
Whole_Picture_Cluster_mode<- function(workdir=WorkingDir(),
cluster_list_file="Cluster.csv",
clustercol="Protein",
Sample_list_file="SMAPLE.csv",
adductslist="M-H",
Threshold=0.15,...
) {
cluster_list<-read.csv(cluster_list_file)
Sample_list<-read.csv(Sample_list)
for (i in 1:length(Sample_list$Sample)){
MALDI_IMAGE_list[as.numeric(Sample_list$Order[i])] <- importImzMl(paste(file.path(Sample_list$Sample[i]),".imzml",sep=""))
}
MPID <-NULL
SMPLIST_D_MODE<-SMPDB[`&`(SMPDB[,"SMPDB.ID"] %in% SMPRESULT[SMPRESULT[,"pScore"]>0.26,"SMPDB.ID"],!(SMPDB[,"HMDB.ID"] %in% unique(SMPLIST[,"HMDB.ID"]))),]
SMPLIST_D_MODE$mass<- t(as.data.frame(lapply(SMPLIST_D_MODE$Formula, getMonomass)))
SMPLIST_D_MODE$moleculeNames<-SMPLIST_D_MODE$Metabolite.Name
for (adduction in 1:length(adductslist)){
SMPLIST_D_MODE$adduct<-adductslist[adduction]
SMPLIST_D_MODE$mz <- SMPLIST_D_MODE$mass+ adducts[adducts[,"Name"]==adductslist[adduction],"Mass"]
if (adduction==1){
META_pathway_ion_image(SMPLIST_D_MODE, SMPRESULT, datafile,MALDI_IMAGE, Image_Type="D-mode", PWpScore=0.26,plotTolerance=0.125,creat_new_file=TRUE)
}else{
META_pathway_ion_image(SMPLIST_D_MODE, SMPRESULT, datafile,MALDI_IMAGE, Image_Type="D-mode", PWpScore=0.26,plotTolerance=0.125,creat_new_file=FALSE)
}
}
options(warn=-1)
SMPLIST_D_MODE_combine<-bind_rows(SMPLIST_D_MODE[,colnames(SMPLIST)],SMPLIST)
options(warn=0)
META_pathway_Whole_Picture(SMPLIST_D_MODE_combine, SMPRESULT, datafile, PWpScore=0.26,pathway_mode="D-mode")
}
#for (j in 1:length(SMPRESULT$SMPDB.ID)){
#if (SMPRESULT$pScore[j]>=PWpScore){
#dir<-paste(datafile,"/",SMPRESULT$SMPDB.ID[j],sep="")
#discovery_moldb_list<-NULL
#if (exists("Map_SMP_MP")==FALSE){Map_SMP_MP <- read.csv("DB/smpdb_pathways.csv")}
#MPID <- sprintf("PW%06d", match(SMPRESULT$SMPDB.ID[j],Map_SMP_MP[,1]))
#PWsvg <- htmlParse(paste("DB/smpdb_svg/",MPID ,".svg",sep=""))
#PWSUBLIST<-SMPLIST[SMPLIST$SMPDB.ID==SMPRESULT$SMPDB.ID[j],c("HMDB.ID","mz","moleculeNames","Metabolite.ID")]
#PW_svg_line<- readLines(paste("DB/smpdb_svg/",MPID ,".svg",sep=""))
#discovery_moldb_cache<-grep('<g id="location" transform=',as.character(PW_svg_line),value=FALSE)
#discovery_moldb_cache<-str_locate(PW_svg_line[discovery_moldb_cache],"data-element-id=")
#for (i in 1 : length(discovery_moldb_cache)){
#temp<-substr(PW_svg_line[discovery_moldb_cache[i]],discovery_moldb_startpos[i,"end"]+2,discovery_moldb_startpos[i,"end"]+2+9)
#discovery_moldb_cache<-grep('data-element-type="compound" data-element-id=',as.character(PW_svg_line),value=FALSE)
#discovery_moldb_startpos<-str_locate(PW_svg_line[discovery_moldb_cache],"data-element-id=")
#discovery_moldb_list<-substr(PW_svg_line[discovery_moldb_cache],discovery_moldb_startpos[,"end"]+2,discovery_moldb_startpos[,"end"]+2+9)
#discovery_moldb_list <- union(discovery_moldb_list,discovery_moldb_list)
#discovery_moldb_list<- fetch_ID_from_DB(discovery_moldb_list,SMPDB,SMPRESULT$SMPDB.ID[j])
##fetch_ID_from SMPDB
fetch_ID_from_DB <- function(Fetch_id_list, SMPDB, SMPDB.ID){
if (SMPDB.ID != "") {for (z in 1: length(Fetch_id_list)){temp3[z,]<-subset(SMPDB , Metabolite.ID == as.character(Fetch_id_list)[z] & SMPDB.ID ==SMPRESULT$SMPDB.ID[j])}}else{
for (z in 1: length(Fetch_id_list)){temp3[z,]<-subset(SMPDB , Metabolite.ID == as.character(Fetch_id_list)[z] )}}
return(temp3)
}
nice_file_create <- function(filename){
ok <- file.create(filename)
if(!ok)
{return(windows_filename2(filename))
}else{return(filename)}
}
Kernel_convolution<-function(Neighbour,resolution=1){
Neighbour<-1+(Neighbour*2*resolution)
if (Neighbour==1){kern = matrix(1, ncol = Neighbour, nrow = Neighbour)}else{
kern = matrix(1/(Neighbour*Neighbour-1), ncol = Neighbour, nrow = Neighbour)
kern[(Neighbour+1)/2,(Neighbour+1)/2]=0}
kern
}
## Basic Function PNG ion image plot
Plot_Ion_image_Png<- function(WKdir, imagefile, png_filename, mz,adduct="M-H", Tolerance= 0.25, title="",Neighbour = 0,Creat_new_file=TRUE,color="black",interpolate =FALSE){
#x11()
#dev.copy(png,paste(dir,"/",png_filename,'.png',sep=""))
Tolerance<-round(Tolerance,digits = 5)
pngfillewrite<-paste(WKdir,"/",png_filename,'.png',sep="")
#pngfillewrite<-nice_file_create(pngfillewrite)
png(paste(WKdir,"/","temp.png",sep=""), bg = "transparent")
try(plotMsiSlice(imagefile ,mz , tolerance=Tolerance, legend=FALSE,colRamp=colorRamp(c("black", "blue", "green", "yellow", "red","#FF00FF","white")),interpolate =interpolate ),silent = TRUE)
dev.off()
try(pngfile<-image_read(paste(WKdir,"/","temp.png",sep="")))
res <- try(pngfile<-image_trim(pngfile),silent = TRUE)
if (class(res) != "try-error"){
#resolution<-imagefile[[1]]@metaData[["imaging"]][["Size"]]["x"]
kern=Kernel_convolution(Neighbour,resolution)
if (title=="D-mode"){
pngfile<-image_border(pngfile, "transparent", "0x30")
pngfile<-image_annotate(pngfile,paste(png_filename,"D-mode","\n",adduct,sprintf("%g",mz),"±",Tolerance),gravity = "south",size = 9)
pngfile<-image_trim(pngfile)
res2 <- try(pngfileoriginal<-image_read(pngfillewrite),silent = TRUE)
if ('&'(class(res2) != "try-error",Creat_new_file==FALSE)){
pngfile<-image_append(c(pngfileoriginal,pngfile))
image_write(pngfile,pngfillewrite)}
if ('&'(class(res2) != "try-error",Creat_new_file==TRUE)){
image_write(pngfile,pngfillewrite)}
if (class(res2) == "try-error"){
image_write(pngfile,pngfillewrite)}
}else{
pngfile<-image_convolve(pngfile, kern)
pngfile<-image_border(pngfile, "transparent", "20x20")
pngfile<-image_annotate(pngfile,paste(png_filename,"\n", adduct,sprintf("%g",mz),"±",Tolerance),gravity = "south",size = 10,color = color)
pngfile<-image_trim(pngfile)
image_write(pngfile,pngfillewrite)
}
}
file.remove(paste(WKdir,"/","temp.png",sep=""))
}
Plot_Ion_image_Png3<- function(WKdir, imagefile, png_filename, mz,adduct="M-H", Tolerance= 0.25, title="",Neighbour = 0,Creat_new_file=TRUE,color="black",interpolate =FALSE){
#x11()
#dev.copy(png,paste(dir,"/",png_filename,'.png',sep=""))
Tolerance<-round(Tolerance,digits = 5)
pngfillewrite<-paste(WKdir,"/",png_filename,'.png',sep="")
#pngfillewrite<-nice_file_create(pngfillewrite)
png(paste(WKdir,"/","temp.png",sep=""), bg = "transparent")
try(plotMsiSlice(imagefile ,mz , tolerance=Tolerance, legend=FALSE,colRamp=colorRamp(c("black", "blue", "green", "yellow", "red","#FF00FF","white")),interpolate =interpolate ),silent = TRUE)
dev.off()
try(pngfile<-image_read(paste(WKdir,"/","temp.png",sep="")))
res <- try(pngfile<-image_trim(pngfile),silent = TRUE)
if (class(res) != "try-error"){
a<-image_attributes(pngfile)
pixelx<-as.numeric(str_split(a[a$property=="png:IHDR.width,height","value"],", ")[[1]][1])
pixely<-as.numeric(str_split(a[a$property=="png:IHDR.width,height","value"],", ")[[1]][2])
try(resolutionx<-pixelx/imagefile[[1]]@metaData[["imaging"]][["size"]]["x"])
try(resolutiony<-pixely/imagefile[[1]]@metaData[["imaging"]][["size"]]["y"])
try(resolution<-sqrt(as.numeric((resolutionx + resolutiony) / 2)))
try(kern<-Kernel_convolution(Neighbour,resolution))
if (title=="D-mode"){
pngfile<-image_border(pngfile, "transparent", "0x30")
pngfile<-image_annotate(pngfile,paste(png_filename,"D-mode","\n",adduct,sprintf("%g",mz),"±",Tolerance),gravity = "south",size = 9)
pngfile<-image_trim(pngfile)
res2 <- try(pngfileoriginal<-image_read(pngfillewrite),silent = TRUE)
if ('&'(class(res2) != "try-error",Creat_new_file==FALSE)){
pngfile<-image_append(c(pngfileoriginal,pngfile))
image_write(pngfile,pngfillewrite)}
if ('&'(class(res2) != "try-error",Creat_new_file==TRUE)){
image_write(pngfile,pngfillewrite)}
if (class(res2) == "try-error"){
image_write(pngfile,pngfillewrite)}
}else{
pngfile<-image_convolve(pngfile, kern)
pngfile<-image_border(pngfile, "transparent", "20x20")
pngfile<-image_annotate(pngfile,paste(png_filename,"\n", adduct,sprintf("%g",mz),"±",Tolerance),gravity = "south",size = 10,color = color)
pngfile<-image_trim(pngfile)
image_write(pngfile,pngfillewrite)
}
}
file.remove(paste(WKdir,"/","temp.png",sep=""))
}
Plot_Ion_image_Png2<- function(WKdir, imagefile, png_filename, mz,adduct="M-H", Tolerance= 0.25, title="",Creat_new_file=TRUE){
#x11()
pngfillewrite<-paste(WKdir,"/",png_filename,'.png',sep="")
#pngfillewrite<-nice_file_create(pngfillewrite)
#dev.copy(png,paste(dir,"/",png_filename,'.png',sep=""))
png(paste(WKdir,"/","temp.png",sep=""))
try(plotMsiSlice(imagefile ,mz , tolerance=Tolerance,colRamp=colorRamp(c("black", "blue", "green", "yellow", "red")), legend=FALSE),silent = TRUE)
dev.off()
try(pngfile<-image_read(paste(WKdir,"/","temp.png",sep="")))
res <- try(pngfile<-image_trim(pngfile),silent = TRUE)
if (class(res) != "try-error"){
if (title=="D-mode"){
pngfile<-image_annotate(pngfile,paste(png_filename,"D-mode","\n",adduct,sprintf("%g",mz),"±",Tolerance),size = 9)
res2 <- try(pngfileoriginal<-image_read(pngfillewrite),silent = TRUE)
if ('&'(class(res2) != "try-error",Creat_new_file==FALSE)){
pngfile<-image_append(c(pngfileoriginal,pngfile))
image_write(pngfile,pngfillewrite)}
if ('&'(class(res2) != "try-error",Creat_new_file==TRUE)){
image_write(pngfile,pngfillewrite)}
if (class(res2) == "try-error"){
image_write(pngfile,pngfillewrite)}
}else{
pngfile<-image_annotate(pngfile,paste(png_filename,"\n", adduct,sprintf("%g",mz),"±",Tolerance),size = 9)
image_write(pngfile,pngfillewrite)
}
}
file.remove(paste(WKdir,"/","temp.png",sep=""))
}
#wrap the getMolecule function
##prepare SMPDB
#prepare_SMPDB<- function(DBpath="DB"){
#SMPDB <- list.files(path=paste(DBpath,"/metapath",sep=""),full.names = TRUE) %>% lapply(read.csv) %>% bind_rows
#write.csv(SMPDB, file=paste(DBpath,"/SMPdb.csv",sep=","))
#if (exists("SMPDB")==FALSE){SMPDB <- read.csv(file="DB/SMPdb.csv", as.is = TRUE)}
#SMPDBFREQ<-as.data.frame(table(SMPDB[,"SMPDB.ID"]))
#write.csv(SMPDBFREQ, file=(paste(DBpath,"/SMPdb_ID_FREQ.csv",sep="")))
#DBpath="DB"
#SMPDB <- read.csv(file=paste(DBpath,"/SMPdb.csv",sep=""), as.is = TRUE)
#SMPDBFREQ <- read.csv(file=paste(DBpath,"/SMPdb_ID_FREQ.csv",sep="") ,col.names= c("","SMPDB.ID","Freq"),as.is = TRUE)
#pathwayDB <- read.csv(file=paste(DBpath,"/smpdb_pathways.csv",sep=""), as.is = TRUE)
#R Python interface
#install.packages("reticulate")
#install.packages("R.matlab")
#}
quiet <- function(x) {
sink(tempfile())
on.exit(sink())
invisible(force(x))
}
Cardinal_old<- function(){
workdir <- WorkingDir()
datafile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.imzML$")
listfile <- list.files(path=workdir,full.names = TRUE,pattern = "\\.csv$")
datafile<-gsub(".imzML", "", datafile)
name <-gsub(paste(workdir,"\\/",sep=""),"",datafile[3])
#imdata <- readImzML(name, folder, as="MSImagingExperiment")
imdata <- readImzML(name, folder, attach.only=TRUE, as="MSImagingExperiment")
summarize(data, sum, .by="pixel")
tmp <- imdata %>%
smoothSignal() %>%
reduceBaseline() %>%
peakPick() %>%
peakFilter() %>%
select(x == 1, y == 1)
process(plot=TRUE,
par=list(layout=c(1,3)),
BPPARAM=SerialParam())
plot(data, pixel=1)
plot(data, coord=list(x=2, y=2))
pattern <- factor(c(0, 0, 2, 2, 0, 0, 0, 0, 0, 0, 2, 2, 0,
0, 0, 0, 0, 0, 0, 1, 2, 2, 0, 0, 0, 0, 0, 2, 1, 1, 2,
2, 0, 0, 0, 0, 0, 1, 2, 2, 2, 2, 0, 0, 0, 0, 1, 2, 2,
2, 2, 2, 0, 0, 0, 0, 2, 2, 2, 2, 2, 2, 2, 0, 0, 0, 2,
2, 0, 0, 0, 0, 0, 0, 2, 2, 0, 0, 0, 0, 0),
levels=c(0,1,2), labels=c("blue", "black", "red"))
tmp <- imdata %>%
smoothSignal() %>%
reduceBaseline() %>%
peakPick() %>%
peakFilter() %>%
select(x == 1, y == 1) %>%
process(plot=TRUE,
par=list(layout=c(1,3)),
BPPARAM=SerialParam())
}
intensity.colors_customize <- function(n = 100, alpha = 1) {
col2 <- rainbow(3*n, alpha=alpha)[(2*n):1]
f <- colorRamp(colors=c("black", "blue", "green", "yellow", "red","#FF00FF","white"))
alpha <- col2rgb(col2, alpha=TRUE)[[4]]
col1 <- sapply(seq(from=0, to=1, length.out=n), function(i) do.call(rgb,
c(as.list(f(i)), maxColorValue=255, alpha=alpha)))
col1
}
Pathway_overview_graphite<-function(){
p_load(graphite,graph )
p_load(org.Hs.eg.db)
p_load(Rgraphviz)
humanReactome <- graphite::pathways("hsapiens", "reactome")
humanmeatbolome <- graphite::pathways("hsapiens", "smpdb")
metab_url <-
url("https://romualdi.bio.unipd.it/wp-uploads/2018/04/Terunuma_metabolite_expr.txt")
mexpr <- read.table(metab_url, header = TRUE, sep = "\t", row.names = NULL,
stringsAsFactors = FALSE)
p <- humanReactome[["Glycolysis"]]
head(nodes(p))
head(edges(p))
head(nodes(p), which = "mixed")
head(edges(p), which = "mixed")
pathwayDatabases()
g <- pathwayGraph(p)
g <- pathwayGraph(p, which = "mixed")
pSymbol <- convertIdentifiers(p, "SYMBOL")
head(nodes(pSymbol))
reactomeSymbol <- convertIdentifiers(humanReactome[1:5], "SYMBOL")
cytoscapePlot(convertIdentifiers(reactome$`Unwinding of DNA`, "symbol"), which = "mixed")
meta_g <- pathwayGraph(p, which = "metabolites")
plot(g)
}
Cardinal_cluster_classification<-function(){
set.seed(1)
skm <- spatialKMeans(imdata, r=2, k=4, method="adaptive")
plot(skm, col=c("pink", "blue", "red","orange","navyblue"), type=c('p','h'), key=FALSE)
image(skm, col=c("pink", "blue", "red","orange","navyblue"), key=FALSE)
set.seed(1)
ssc <- spatialShrunkenCentroids(imdata, r=2, k=5, s=3, method="gaussian")
plot(ssc, col=c("pink", "blue", "red","orange","navyblue"), type=c('p','h'), key=FALSE)
image(ssc, col=c("pink", "blue", "red","orange","navyblue"), key=FALSE)
}
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.