Mercurial > repos > marie-tremblay-metatoul > nmr_alignment
comparison nmr_alignement/NmrAlignment_wrapper.R @ 3:f3ec6799c435 draft
Uploaded
author | marie-tremblay-metatoul |
---|---|
date | Wed, 12 Apr 2017 04:32:23 -0400 |
parents | 908e1345d7ca |
children | 265ee538e120 |
comparison
equal
deleted
inserted
replaced
2:908e1345d7ca | 3:f3ec6799c435 |
---|---|
15 | 15 |
16 ##------------------------------ | 16 ##------------------------------ |
17 ## Libraries loading | 17 ## Libraries loading |
18 ##------------------------------ | 18 ##------------------------------ |
19 # ParseCommandArgs function | 19 # ParseCommandArgs function |
20 library(batch) | 20 library(batch) |
21 # Alignment | 21 # Alignment |
22 library(speaq) | 22 library(speaq) |
23 | 23 |
24 | 24 |
25 # R script call | 25 # R script call |
49 ## Script | 49 ## Script |
50 ##------------------------------ | 50 ##------------------------------ |
51 if(!runExampleL) | 51 if(!runExampleL) |
52 argLs <- parseCommandArgs(evaluate=FALSE) | 52 argLs <- parseCommandArgs(evaluate=FALSE) |
53 | 53 |
54 | 54 |
55 ## Parameters Loading | 55 ## Parameters Loading |
56 ##------------------- | 56 ##------------------- |
57 # Inputs | 57 # Inputs |
58 ## Library of spectra to align | 58 ## Library of spectra to align |
59 if (!is.null(argLs[["zipfile"]])){ | 59 if (!is.null(argLs[["zipfile"]])){ |
60 fileType="zip" | 60 fileType="zip" |
61 zipfile= argLs[["zipfile"]] | 61 zipfile= argLs[["zipfile"]] |
62 directory=unzip(zipfile, list=F) | 62 directory=unzip(zipfile, list=F) |
63 directory=paste(getwd(),strsplit(directory[1],"/")[[1]][2],sep="/") | 63 directory=paste(getwd(),strsplit(directory[1],"/")[[1]][2],sep="/") |
64 } else if (!is.null(argLs[["library"]])){ | |
65 fileType="zip" | |
66 directory=argLs[["library"]] | |
67 if(!file.exists(directory)){ | |
68 error_message=paste("Cannot access the directory :",directory,".Please verify if the directory exists or not.") | |
69 print(error_message) | |
70 stop(error_message) | |
71 } | |
72 } else if (!is.null(argLs[["tsvfile"]])){ | 64 } else if (!is.null(argLs[["tsvfile"]])){ |
73 fileType="tsv" | 65 fileType="tsv" |
74 directory <- read.table(argLs[["tsvfile"]],check.names=FALSE,header=TRUE,sep="\t") | 66 directory <- read.table(argLs[["tsvfile"]],check.names=FALSE,header=TRUE,sep="\t") |
75 } | 67 } |
76 | 68 |
77 | 69 |
78 ## Spectral width | 70 ## Spectral width |
79 leftBorder <- argLs[["left_border"]] | 71 leftBorder <- argLs[["left_border"]] |
80 rightBorder <- argLs[["right_border"]] | 72 rightBorder <- argLs[["right_border"]] |
81 | 73 |
82 ##Exclusion zone(s) | 74 ##Exclusion zone(s) |
83 exclusionZones <- argLs[["zone_exclusion_choices.choice"]] | 75 exclusionZones <- argLs[["zone_exclusion_choices.choice"]] |
84 exclusionZonesBorders <- NULL | 76 exclusionZonesBorders <- NULL |
85 if (!is.null(argLs$zone_exclusion_left)) | 77 if (!is.null(argLs$zone_exclusion_left)) |
86 { | 78 { |
87 for(i in which(names(argLs)=="zone_exclusion_left")) | 79 for(i in which(names(argLs)=="zone_exclusion_left")) |
88 { | 80 { |
89 exclusionZonesBorders <- c(exclusionZonesBorders,list(c(argLs[[i]],argLs[[i+1]]))) | 81 exclusionZonesBorders <- c(exclusionZonesBorders,list(c(argLs[[i]],argLs[[i+1]]))) |
90 } | 82 } |
91 } | 83 } |
92 | 84 |
93 ## Reference spectrum | 85 ## Reference spectrum |
94 reference <- argLs[["reference"]] | 86 reference <- argLs[["reference"]] |
95 | 87 |
96 ## Size of a small nDivRange | 88 ## Size of a small nDivRange |
97 nDivRange <- argLs[["nDivRange"]] | 89 nDivRange <- argLs[["nDivRange"]] |
98 | 90 |
99 ## Intensity threshold for peak removal | 91 ## Intensity threshold for peak removal |
100 baselineThresh <- argLs[["baselineThresh"]] | 92 baselineThresh <- argLs[["baselineThresh"]] |
101 | 93 |
102 | 94 |
103 # Outputs | 95 # Outputs |
109 ## Checking arguments | 101 ## Checking arguments |
110 ##------------------- | 102 ##------------------- |
111 error.stock <- "\n" | 103 error.stock <- "\n" |
112 if(length(error.stock) > 1) | 104 if(length(error.stock) > 1) |
113 stop(error.stock) | 105 stop(error.stock) |
114 | 106 |
115 | 107 |
116 ## Computation | 108 ## Computation |
117 ##------------ | 109 ##------------ |
118 directory.alignement <- nmr.alignment(fileType=fileType,directory=directory,leftBorder=leftBorder,rightBorder=rightBorder,exclusionZones=exclusionZones, | 110 directory.alignement <- nmr.alignment(fileType=fileType,directory=directory,leftBorder=leftBorder,rightBorder=rightBorder,exclusionZones=exclusionZones, |
119 exclusionZonesBorders=exclusionZonesBorders, reference=reference, nDivRange=nDivRange, | 111 exclusionZonesBorders=exclusionZonesBorders, reference=reference, nDivRange=nDivRange, |
120 baselineThresh=baselineThresh, maxshift=50, verbose=FALSE) | 112 baselineThresh=baselineThresh, maxshift=50, verbose=FALSE) |
121 directory.raw <- directory.alignement[[1]] | 113 directory.raw <- directory.alignement[[1]] |
122 directory.aligned <- directory.alignement[[2]] | 114 directory.aligned <- directory.alignement[[2]] |
123 | 115 |
124 ## Saving | 116 ## Saving |
125 ##------- | 117 ##------- |
126 # Aligned spectra | 118 # Aligned spectra |
127 t.directory.aligned <- t(directory.aligned) | 119 t.directory.aligned <- t(directory.aligned) |
128 rownames(t.directory.aligned) <- colnames(directory.aligned) | 120 rownames(t.directory.aligned) <- colnames(directory.aligned) |
129 # colnames(t.directory.aligned) <- c("Bucket",colnames(t.directory.aligned)) | 121 # colnames(t.directory.aligned) <- c("Bucket",colnames(t.directory.aligned)) |
130 write.table(t.directory.aligned,file=alignedSpectra,row.names=TRUE,quote=FALSE,sep="\t") | 122 write.table(t.directory.aligned,file=alignedSpectra,row.names=TRUE,quote=FALSE,sep="\t") |
131 | 123 |
132 | 124 |
133 excludedZone <- NULL | 125 excludedZone <- NULL |
134 for (c in 1:length(exclusionZonesBorders)) | 126 for (c in 1:length(exclusionZonesBorders)) |
135 { | 127 { |
162 { | 154 { |
163 drawSpec(raw.spectra[,(which(round(as.numeric(colnames(raw.spectra)),2) == excludedZone[n])[1]):(which(round(as.numeric(colnames(raw.spectra)),2) == excludedZone[n-1])[1])],xlab="", ylab="Raw spectra", main="") | 155 drawSpec(raw.spectra[,(which(round(as.numeric(colnames(raw.spectra)),2) == excludedZone[n])[1]):(which(round(as.numeric(colnames(raw.spectra)),2) == excludedZone[n-1])[1])],xlab="", ylab="Raw spectra", main="") |
164 drawSpec(aligned.spectra[,(which(round(as.numeric(colnames(aligned.spectra)),2) == excludedZone[n])[1]):(which(round(as.numeric(colnames(aligned.spectra)),2) == excludedZone[n-1])[1])],xlab="", ylab="Aligned spectra", main="") | 156 drawSpec(aligned.spectra[,(which(round(as.numeric(colnames(aligned.spectra)),2) == excludedZone[n])[1]):(which(round(as.numeric(colnames(aligned.spectra)),2) == excludedZone[n-1])[1])],xlab="", ylab="Aligned spectra", main="") |
165 n <- n - 2 | 157 n <- n - 2 |
166 } | 158 } |
167 | 159 |
168 drawSpec(raw.spectra[,(which(round(as.numeric(colnames(raw.spectra)),2) == excludedZone[1])[1]):ncol(raw.spectra)],xlab="", ylab="Raw spectra", main="") | 160 drawSpec(raw.spectra[,(which(round(as.numeric(colnames(raw.spectra)),2) == excludedZone[1])[1]):ncol(raw.spectra)],xlab="", ylab="Raw spectra", main="") |
169 drawSpec(aligned.spectra[,(which(round(as.numeric(colnames(aligned.spectra)),2) == excludedZone[1])[1]):ncol(aligned.spectra)],xlab="", ylab="Aligned spectra", main="") | 161 drawSpec(aligned.spectra[,(which(round(as.numeric(colnames(aligned.spectra)),2) == excludedZone[1])[1]):ncol(aligned.spectra)],xlab="", ylab="Aligned spectra", main="") |
170 } | 162 } |
171 drawSpec(raw.spectra[,(which(round(as.numeric(colnames(raw.spectra)),2) == 2.4)[1]):(which(round(as.numeric(colnames(raw.spectra)),2) == 2.8)[1])],xlab="", ylab="Raw spectra", main="") | 163 drawSpec(raw.spectra[,(which(round(as.numeric(colnames(raw.spectra)),2) == 2.4)[1]):(which(round(as.numeric(colnames(raw.spectra)),2) == 2.8)[1])],xlab="", ylab="Raw spectra", main="") |
172 drawSpec(aligned.spectra[,(which(round(as.numeric(colnames(aligned.spectra)),2) == 2.4)[1]):(which(round(as.numeric(colnames(aligned.spectra)),2) == 2.8)[1])],xlab="", ylab="Aligned spectra", main="") | 164 drawSpec(aligned.spectra[,(which(round(as.numeric(colnames(aligned.spectra)),2) == 2.4)[1]):(which(round(as.numeric(colnames(aligned.spectra)),2) == 2.8)[1])],xlab="", ylab="Aligned spectra", main="") |