annotate blupcal_wrapper.R @ 0:45d215f2be74 draft

Uploaded
author dereeper
date Sat, 29 Dec 2018 18:44:05 -0500
parents
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
0
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
1 #!/usr/bin/Rscript
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
2
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
3 options(show.error.messages = F,
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
4 error = function() {
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
5 cat(geterrmessage(), file = stderr())
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
6 q("no", 1, F) } )
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
7 loc <- Sys.setlocale("LC_MESSAGES", "en_US.UTF-8")
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
8 library("getopt")
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
9
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
10 options(stringAsfactors = FALSE, useFancyQuotes = FALSE)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
11 options(warn = -1)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
12 args <- commandArgs(trailingOnly = TRUE)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
13 option_specification <- matrix(c(
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
14 "tool_directory", "a", 1, "character"
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
15 ,"tabular_file", "b", 1, "character"
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
16 ,"replication_vector_column_index", "c", 1, "integer"
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
17 ,"genotype_vector_column_index", "d", 1, "integer"
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
18 ,"y_vector_column_index", "e", 1, "integer"
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
19 ,"design", "f", 1, "character"
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
20 ,"block_vector_column_index", "g", 1, "integer"
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
21 ,"summarize_by", "h", 1, "character"
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
22 ,"summarize_by_vector_column_index", "i", 1, "integer"
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
23 ,"group_variable_1", "j", 1, "character"
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
24 ,"group_variable_1_vector_column_index", "k", 1, "integer"
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
25 ,"group_variable_2", "l", 1, "character"
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
26 ,"group_variable_2_vector_column_index", "m", 1, "integer"
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
27 ,"output_file_path", "n", 1, "character"
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
28 ,"png_plots_file_path", "o", 1, "character"
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
29 ,"png_histograms_file_path", "p", 1, "character"
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
30 ,"pdf_plots_file_path", "q", 1, "character"
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
31 ,"pdf_histograms_file_path", "r", 1, "character"
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
32 ), byrow = TRUE, ncol = 4)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
33 options <- getopt(option_specification)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
34
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
35 cat("\nTool Directory: ", options$tool_directory)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
36 cat("\nTabular File: ", options$tabular_file)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
37 cat("\nReplication Vector Column Index: ", options$replication_vector_column_index)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
38 cat("\nGenotype Vector Column Index: ", options$genotype_vector_column_index)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
39 cat("\nY Vector Column Index: ", options$y_vector_column_index)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
40 cat("\nDesign: ", options$design)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
41 cat("\nBlock Vector Column Index: ", options$block_vector_column_index)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
42 cat("\nSummarize By: ", options$summarize_by)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
43 cat("\nSummarize By Vector Column Index: ", options$summarize_by_vector_column_index)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
44 cat("\nGroup Variable #1: ", options$group_variable_1)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
45 cat("\nGroup Variable #1 Vector Column Index: ", options$group_variable_1_vector_column_index)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
46 cat("\nGroup Variable #2: ", options$group_variable_2)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
47 cat("\nGroup Variable #2 Vector Column Index: ", options$group_variable_2_vector_column_index)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
48 cat("\nOutput file path: ", options$output_file_path)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
49 cat("\nPNG plots file path: ", options$png_plots_file_path)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
50 cat("\nPNG histograms file path: ", options$png_histograms_file_path)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
51 cat("\nPDF plots file path: ", options$pdf_plots_file_path)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
52 cat("\nPDF histograms file path: ", options$pdf_histograms_file_path)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
53 cat("\n\n")
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
54
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
55 tabular_data <- read.table(file = options$tabular_file,
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
56 header = TRUE,
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
57 sep = "\t",
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
58 stringsAsFactors = FALSE,
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
59 strip.white = TRUE,
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
60 na.strings = ".")
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
61 #tabular_data
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
62
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
63 column_names <- colnames(tabular_data)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
64 cat("Column names:\n")
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
65 column_names
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
66 cat("\n\n")
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
67
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
68 # Replication
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
69 replication_vector_header <- column_names[options$replication_vector_column_index]
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
70 cat("\nreplication vector header: ", replication_vector_header)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
71
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
72 # Genotype
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
73 genotype_vector_header <- column_names[options$genotype_vector_column_index]
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
74 cat("\ngenotype vector header: ", genotype_vector_header)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
75
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
76 # Y
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
77 y_vector_header <- column_names[options$y_vector_column_index]
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
78 cat("\ny vector header: ", y_vector_header)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
79
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
80 # Design
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
81 design <- options$design
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
82 cat("\ndesign: ", design)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
83
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
84 # Block
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
85 if (design == "rcbd") {
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
86 block_vector_header <- NULL
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
87 } else {
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
88 block_vector_header <- column_names[options$block_vector_column_index]
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
89 }
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
90 cat("\nblock vector header: ", block_vector_header)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
91
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
92 # Summarize By
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
93 if (options$summarize_by == "false") {
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
94 summarize_by_vector_header <- NULL
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
95 } else {
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
96 summarize_by_vector_header <- column_names[options$summarize_by_vector_column_index]
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
97 }
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
98 cat("\nsummarize by vector header: ", summarize_by_vector_header)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
99
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
100 # Group Variable #1
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
101 if (options$group_variable_1 == "false") {
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
102 group_variable_1_vector_header <- NULL
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
103 } else {
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
104 group_variable_1_vector_header <- column_names[options$group_variable_1_vector_column_index]
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
105 }
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
106 cat("\ngroup variable 1 vector header: ", group_variable_1_vector_header)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
107
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
108 # Group Variable #2
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
109 if (options$group_variable_2 == "false") {
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
110 group_variable_2_vector_header <- NULL
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
111 } else {
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
112 group_variable_2_vector_header <- column_names[options$group_variable_2_vector_column_index]
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
113 }
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
114 cat("\ngroup variable 2 vector header: ", group_variable_2_vector_header)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
115
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
116 cat("\n\n")
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
117
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
118 source(paste(options$tool_directory, "/blupcal.R", sep = ""))
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
119
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
120 fit <- blupcal(data = tabular_data,
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
121 Replication = replication_vector_header,
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
122 Genotype = genotype_vector_header,
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
123 y = y_vector_header,
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
124 design = design,
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
125 Block = block_vector_header,
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
126 summarizeby = summarize_by_vector_header,
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
127 groupvar1 = group_variable_1_vector_header,
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
128 groupvar2 = group_variable_2_vector_header)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
129
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
130 cat("\n\n")
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
131
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
132 output_file_path <- options$output_file_path
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
133
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
134 cat("\nfit: ")
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
135 cat("\n")
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
136 fit
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
137 cat("\n")
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
138
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
139 #fit$gid <- as.numeric(levels(fit$gid))[fit$gid]
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
140 if (options$summarize_by == "false") {
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
141 fit_view <- fit
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
142 } else {
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
143 fit_view <- fit[c(1, 7, 2, 3, 4, 5, 6)]
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
144 }
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
145
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
146 formatted_table <- lapply(fit_view, function(.col) {
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
147 if (is.numeric(.col)) {
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
148 return(as.numeric(sprintf("%.3f", .col)))
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
149 } else {
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
150 return(.col)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
151 }
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
152 })
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
153
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
154 write.table(formatted_table,
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
155 file = output_file_path,
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
156 sep = "\t",
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
157 row.names = FALSE,
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
158 quote = FALSE)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
159
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
160 png_plots_file_path <-options$png_plots_file_path
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
161 if (!(is.null(png_plots_file_path))) {
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
162 png(png_plots_file_path)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
163 plot(fit)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
164 dev.off()
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
165 }
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
166
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
167 png_histograms_file_path <-options$png_histograms_file_path
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
168 if (!(is.null(png_histograms_file_path))) {
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
169 png(png_histograms_file_path)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
170 plot2.blupcal(fit)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
171 dev.off()
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
172 }
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
173
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
174 pdf_plots_file_path <-options$pdf_plots_file_path
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
175 if (!(is.null(pdf_plots_file_path))) {
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
176 pdf(pdf_plots_file_path)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
177 plot(fit)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
178 dev.off()
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
179 }
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
180
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
181 pdf_histograms_file_path <-options$pdf_histograms_file_path
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
182 if (!(is.null(pdf_histograms_file_path))) {
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
183 pdf(pdf_histograms_file_path)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
184 plot2.blupcal(fit)
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
185 dev.off()
45d215f2be74 Uploaded
dereeper
parents:
diff changeset
186 }