Mercurial > repos > george-weingart > maaslin

comparison src/lib/scriptBiplotTSV.R @ 0:e0b5980139d9

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
maaslin
author george-weingart
date Tue, 13 May 2014 22:00:40 -0400
parents
children
comparison
equal deleted inserted replaced
-1:000000000000 0:e0b5980139d9
1 #!/usr/bin/env Rscript
2
3 library(vegan)
4 library(optparse)
5
6 funcGetCentroidForMetadatum <- function(
7 ### Given a binary metadatum, calculate the centroid of the samples associated with the metadata value of 1
8 # 1. Get all samples that have the metadata value of 1
9 # 2. Get the x and y coordinates of the selected samples
10 # 3. Get the median value for the x and ys
11 # 4. Return those coordinates as the centroid's X and Y value
12 vfMetadata,
13 ### Logical or integer (0,1) vector, TRUE or 1 values indicate correspoinding samples in the
14 ### mSamplePoints which will be used to define the centroid
15 mSamplePoints
16 ### Coordinates (columns;n=2) of samples (rows) corresponding to the vfMetadata
17 ){
18 # Check the lengths which should be equal
19 if(length(vfMetadata)!=nrow(mSamplePoints))
20 {
21 print(paste("funcGetCentroidForMetadata::Error: Should have received metadata and samples of the same length, received metadata length ",length(vfMetadata)," and sample ",nrow(mSamplePoints)," length.",sep=""))
22 return( FALSE )
23 }
24
25 # Get all the samples that have the metadata value of 1
26 viMetadataSamples = which(as.integer(vfMetadata)==1)
27
28 # Get the x and y coordinates for the selected samples
29 mSelectedPoints = mSamplePoints[viMetadataSamples,]
30
31 # Get the median value for the x and the ys
32 if(!is.null(nrow(mSelectedPoints)))
33 {
34 return( list(x=median(mSelectedPoints[,1],na.rm = TRUE),y=median(mSelectedPoints[,2],na.rm = TRUE)) )
35 } else {
36 return( list(x=mSelectedPoints[1],y=mSelectedPoints[2]) )
37 }
38 }
39
40 funcGetMaximumForMetadatum <- function(
41 ### Given a continuous metadata
42 ### 1. Use the x and ys from mSamplePoints for coordinates and the metadata value as a height (z)
43 ### 2. Use lowess to smooth the landscape
44 ### 3. Take the maximum of the landscape
45 ### 4. Return the coordiantes for the maximum as the centroid
46 vdMetadata,
47 ### Continuous (numeric or integer) metadata
48 mSamplePoints
49 ### Coordinates (columns;n=2) of samples (rows) corresponding to the vfMetadata
50 ){
51 # Work with data frame
52 if(class(mSamplePoints)=="matrix")
53 {
54 mSamplePoints = data.frame(mSamplePoints)
55 }
56 # Check the lengths of the dataframes and the metadata
57 if(length(vdMetadata)!=nrow(mSamplePoints))
58 {
59 print(paste("funcGetMaximumForMetadatum::Error: Should have received metadata and samples of the same length, received metadata length ",length(vdMetadata)," and sample ",nrow(mSamplePoints)," length.",sep=""))
60 return( FALSE )
61 }
62
63 # Add the metadata value to the points
64 mSamplePoints[3] = vdMetadata
65 names(mSamplePoints) = c("x","y","z")
66
67 # Create lowess to smooth the surface
68 # And calculate the fitted heights
69 # x = sample coordinate 1
70 # y = sample coordinate 2
71 # z = metadata value
72 loessSamples = loess(z~x*y, data=mSamplePoints, degree = 1, normalize = FALSE, na.action=na.omit)
73
74 # Naively get the max
75 vdCoordinates = loessSamples$x[which(loessSamples$y==max(loessSamples$y)),]
76 return(list(lsmod = loessSamples, x=vdCoordinates[1],y=vdCoordinates[2]))
77 }
78
79 funcMakeShapes <- function(
80 ### Takes care of defining shapes for the plot
81 dfInput,
82 ### Data frame of metadata measurements
83 sShapeBy,
84 ### The metadata to shape by
85 sShapes,
86 ### List of custom metadata (per level if factor).
87 ### Should correspond to the number of levels in shapeBy; the format is level:shape,level:shape for example HighLuminosity:14,LowLuminosity:2,HighPH:10,LowPH:18
88 cDefaultShape
89 ### Shape to default to if custom shapes are not used
90 ){
91 lShapes = list()
92 vsShapeValues = c()
93 vsShapeShapes = c()
94 vsShapes = c()
95 sMetadataId = sShapeBy
96
97 # Set default shape, color, and color ranges
98 if(!is.null(cDefaultShape))
99 {
100 # Default shape should be an int for the int pch options
101 if(!is.na(as.integer(cDefaultShape)))
102 {
103 cDefaultShape = as.integer(cDefaultShape)
104 }
105 } else {
106 cDefaultShape = 16
107 }
108
109 # Make shapes
110 vsShapes = rep(cDefaultShape,nrow(dfInput))
111
112 if(!is.null(sMetadataId))
113 {
114 if(is.null(sShapes))
115 {
116 vsShapeValues = unique(dfInput[[sMetadataId]])
117 vsShapeShapes = 1:length(vsShapeValues)
118 } else {
119 # Put the markers in the order of the values)
120 vsShapeBy = unlist(strsplit(sShapes,","))
121 for(sShapeBy in vsShapeBy)
122 {
123 vsShapeByPieces = unlist(strsplit(sShapeBy,":"))
124 lShapes[vsShapeByPieces[1]] = as.integer(vsShapeByPieces[2])
125 }
126 vsShapeValues = names(lShapes)
127 }
128
129 # Shapes in the correct order
130 if(!is.null(sShapes))
131 {
132 vsShapeShapes = unlist(lapply(vsShapeValues,function(x) lShapes[[x]]))
133 }
134 vsShapeValues = paste(vsShapeValues)
135
136 # Make the list of shapes
137 for(iShape in 1:length(vsShapeValues))
138 {
139 vsShapes[which(paste(dfInput[[sMetadataId]])==vsShapeValues[iShape])]=vsShapeShapes[iShape]
140 }
141
142 # If they are all numeric characters, make numeric
143 viIntNas = which(is.na(as.integer(vsShapes)))
144 viNas = which(is.na(vsShapes))
145 if(length(setdiff(viIntNas,viNas))==0)
146 {
147 vsShapes = as.integer(vsShapes)
148 } else {
149 print("funcMakeShapes::Error: Please supply numbers 1-25 for shape in the -y,--shapeBy option")
150 vsShapeValues = c()
151 vsShapeShapes = c()
152 }
153 }
154 return(list(PlotShapes=vsShapes,Values=vsShapeValues,Shapes=vsShapeShapes,ID=sMetadataId,DefaultShape=cDefaultShape))
155 }
156
157 ### Global defaults
158 c_sDefaultColorBy = NULL
159 c_sDefaultColorRange = "orange,cyan"
160 c_sDefaultTextColor = "black"
161 c_sDefaultArrowColor = "cyan"
162 c_sDefaultArrowTextColor = "Blue"
163 c_sDefaultNAColor = "grey"
164 c_sDefaultShapeBy = NULL
165 c_sDefaultShapes = NULL
166 c_sDefaultMarker = "16"
167 c_sDefaultRotateByMetadata = NULL
168 c_sDefaultResizeArrow = 1
169 c_sDefaultTitle = "Custom Biplot of Bugs and Samples - Metadata Plotted with Centroids"
170 c_sDefaultOutputFile = NULL
171
172 ### Create command line argument parser
173 pArgs <- OptionParser( usage = "%prog last_metadata input.tsv" )
174
175 # Selecting features to plot
176 pArgs <- add_option( pArgs, c("-b", "--bugs"), type="character", action="store", default=NULL, dest="sBugs", metavar="BugsToPlot", help="Comma delimited list of data to plot as text. Bug|1,Bug|2")
177 pArgs <- add_option( pArgs, c("-m", "--metadata"), type="character", action="store", default=NULL, dest="sMetadata", metavar="MetadataToPlot", help="Comma delimited list of metadata to plot as arrows. metadata1,metadata2,metadata3")
178
179 # Colors
180 pArgs <- add_option( pArgs, c("-c", "--colorBy"), type="character", action="store", default=c_sDefaultColorBy, dest="sColorBy", metavar="MetadataToColorBy", help="The id of the metadatum to use to make the marker colors. Expected to be a continuous metadata.")
181 pArgs <- add_option( pArgs, c("-r", "--colorRange"), type="character", action="store", default=c_sDefaultColorRange, dest="sColorRange", metavar="ColorRange", help=paste("Colors used to color the samples; a gradient will be formed between the color.Default=", c_sDefaultColorRange))
182 pArgs <- add_option( pArgs, c("-t", "--textColor"), type="character", action="store", default=c_sDefaultTextColor, dest="sTextColor", metavar="TextColor", help=paste("The color bug features will be plotted with as text. Default =", c_sDefaultTextColor))
183 pArgs <- add_option( pArgs, c("-a", "--arrowColor"), type="character", action="store", default=c_sDefaultArrowColor, dest="sArrowColor", metavar="ArrowColor", help=paste("The color metadata features will be plotted with as an arrow and text. Default", c_sDefaultArrowColor))
184 pArgs <- add_option( pArgs, c("-w", "--arrowTextColor"), type="character", action="store", default=c_sDefaultArrowTextColor, dest="sArrowTextColor", metavar="ArrowTextColor", help=paste("The color for the metadata text ploted by the head of the metadata arrow. Default", c_sDefaultArrowTextColor))
185 pArgs <- add_option(pArgs, c("-n","--plotNAColor"), type="character", action="store", default=c_sDefaultNAColor, dest="sPlotNAColor", metavar="PlotNAColor", help=paste("Plot NA values as this color. Example -n", c_sDefaultNAColor))
186
187 # Shapes
188 pArgs <- add_option( pArgs, c("-y", "--shapeby"), type="character", action="store", default=c_sDefaultShapeBy, dest="sShapeBy", metavar="MetadataToShapeBy", help="The metadata to use to make marker shapes. Expected to be a discrete metadatum. An example would be -y Environment")
189 pArgs <- add_option( pArgs, c("-s", "--shapes"), type="character", action="store", default=c_sDefaultShapes, dest="sShapes", metavar="ShapesForPlotting", help="This is to be used to specify the shapes to use for plotting. Can use numbers recognized by R as shapes (see pch). Should correspond to the number of levels in shapeBy; the format is level:shape,level:shape for example HighLuminosity:14,LowLuminosity:2,HighPH:10,LowPH:18 . Need to specify -y/--shapeBy for this option to work.")
190 pArgs <- add_option( pArgs, c("-d", "--defaultMarker"), type="character", action="store", default=c_sDefaultMarker, dest="sDefaultMarker", metavar="DefaultColorMarker", help="Default shape for markers which are not otherwise indicated in --shapes, can be used for unspecified values or NA. Must not be a shape in --shapes.")
191
192 # Plot manipulations
193 pArgs <- add_option( pArgs, c("-e","--rotateByMetadata"), type="character", action="store", default=c_sDefaultRotateByMetadata, dest="sRotateByMetadata", metavar="RotateByMetadata", help="Rotate the ordination by a metadata. Give both the metadata and value to weight it by. The larger the weight, the more the ordination is influenced by the metadata. If the metadata is continuous, use the metadata id; if the metadata is discrete, the ordination will be by one of the levels so use the metadata ID and level seperated by a '_'. Discrete example -e Environment_HighLumninosity,100 ; Continuous example -e Environment,100 .")
194 pArgs <- add_option( pArgs, c("-z","--resizeArrow"), type="numeric", action="store", default=c_sDefaultResizeArrow, dest="dResizeArrow", metavar="ArrowScaleFactor", help="A constant to multiple the length of the arrow to expand or shorten all arrows together. This will not change the angle of the arrow nor the relative length of arrows to each other.")
195
196 # Misc
197 pArgs <- add_option( pArgs, c("-i", "--title"), type="character", action="store", default=c_sDefaultTitle, dest="sTitle", metavar="Title", help="This is the title text to add to the plot.")
198 pArgs <- add_option( pArgs, c("-o", "--outputfile"), type="character", action="store", default=c_sDefaultOutputFile, dest="sOutputFileName", metavar="OutputFile", help="This is the name for the output pdf file. If an output file is not given, an output file name is made based on the input file name.")
199
200 funcDoBiplot <- function(
201 ### Perform biplot. Samples are markers, bugs are text, and metadata are text with arrows. Markers and bugs are dtermined usiing NMDS and Bray-Curtis dissimilarity. Metadata are placed on the ordination in one of two ways: 1. Factor data - for each level take the ordination points for the samples that have that level and plot the metadata text at the average orindation point. 2. For continuous data - make a landscape (x and y form ordination of the points) and z (height) as the metadata value. Use a lowess line to get the fitted values for z and take the max of the landscape. Plot the metadata text at that smoothed max.
202 sBugs,
203 ### Comma delimited list of data to plot as text. Bug|1,Bug|2
204 sMetadata,
205 ### Comma delimited list of metadata to plot as arrows. metadata1,metadata2,metadata3.
206 sColorBy = c_sDefaultColorBy,
207 ### The id of the metadatum to use to make the marker colors. Expected to be a continuous metadata.
208 sColorRange = c_sDefaultColorRange,
209 ### Colors used to color the samples; a gradient will be formed between the color. Example orange,cyan
210 sTextColor = c_sDefaultTextColor,
211 ### The color bug features will be plotted with as text. Example black
212 sArrowColor = c_sDefaultArrowColor,
213 ### The color metadata features will be plotted with as an arrow and text. Example cyan
214 sArrowTextColor = c_sDefaultArrowTextColor,
215 ### The color for the metadata text ploted by the head of the metadata arrow. Example Blue
216 sPlotNAColor = c_sDefaultNAColor,
217 ### Plot NA values as this color. Example grey
218 sShapeBy = c_sDefaultShapeBy,
219 ### The metadata to use to make marker shapes. Expected to be a discrete metadatum.
220 sShapes = c_sDefaultShapes,
221 ### This is to be used to specify the shapes to use for plotting. Can use numbers recognized by R as shapes (see pch). Should correspond to the number of levels in shapeBy; the format is level:shape,level:shape for example HighLuminosity:14,LowLuminosity:2,HighPH:10,LowPH:18 . Works with sShapesBy.
222 sDefaultMarker = c_sDefaultMarker,
223 ### The default marker shape to use if shapes are not otherwise indicated.
224 sRotateByMetadata = c_sDefaultRotateByMetadata,
225 ### Metadata and value to rotate by. example Environment_HighLumninosity,100
226 dResizeArrow = c_sDefaultResizeArrow,
227 ### Scale factor to resize tthe metadata arrows
228 sTitle = c_sDefaultTitle,
229 ### The title for the figure.
230 sInputFileName,
231 ### File to input (tsv file: tab separated, row = sample file)
232 sLastMetadata,
233 ### Last metadata that seperates data and metadata
234 sOutputFileName = c_sDefaultOutputFile
235 ### The file name to save the figure.
236 ){
237 print("IN Biplot")
238 # Define the colors
239 vsColorRange = c("blue","orange")
240 cDefaultColor = "black"
241 if(!is.null(sColorRange))
242 {
243 vsColorRange = unlist(strsplit(sColorRange,","))
244 }
245
246 # List of bugs to plot
247 # If there is a list it needs to be more than one.
248 vsBugsToPlot = c()
249 if(!is.null(sBugs))
250 {
251 vsBugsToPlot = unlist(strsplit(sBugs,","))
252 }
253
254 print("vsBugsToPlot")
255 print(vsBugsToPlot)
256 # Metadata to plot
257 vsMetadata = c()
258 if(!is.null(sMetadata))
259 {
260 vsMetadata = unlist(strsplit(sMetadata,","))
261 }
262
263 print("vsMetadata")
264 print(vsMetadata)
265 ### Load table
266 if(class(sInputFileName)=="character")
267 {
268 dfInput = read.table(sInputFileName, sep = "\t", header=TRUE)
269 names(dfInput) = unlist(lapply(names(dfInput),function(x) gsub(".","|",x,fixed=TRUE)))
270 row.names(dfInput) = dfInput[,1]
271 dfInput = dfInput[-1]
272 } else {dfInput = sInputFileName}
273
274 ### Get positions of all metadata or all data
275 iLastMetadata = which(names(dfInput)==sLastMetadata)
276 viMetadata = 1:iLastMetadata
277 viData = (iLastMetadata+1):ncol(dfInput)
278
279 ### Dummy the metadata if discontinuous
280 ### Leave the continous metadata alone but include
281 listMetadata = list()
282 vsRowNames = c()
283 viContinuousMetadata = c()
284 for(i in viMetadata)
285 {
286 print( names( dfInput )[i] )
287 vCurMetadata = unlist(dfInput[i])
288 if( ( is.numeric(vCurMetadata)||is.integer(vCurMetadata) ) && ( length( unique( vCurMetadata ) ) >= c_iNonFactorLevelThreshold ) )
289 {
290 vCurMetadata[which(is.na(vCurMetadata))] = mean(vCurMetadata,na.rm=TRUE)
291 listMetadata[[length(listMetadata)+1]] = vCurMetadata
292 vsRowNames = c(vsRowNames,names(dfInput)[i])
293 viContinuousMetadata = c(viContinuousMetadata,length(listMetadata))
294 } else {
295 vCurMetadata = as.factor(vCurMetadata)
296 vsLevels = levels(vCurMetadata)
297 for(sLevel in vsLevels)
298 {
299 vNewMetadata = rep(0,length(vCurMetadata))
300 vNewMetadata[which(vCurMetadata == sLevel)] = 1
301 listMetadata[[length(listMetadata)+1]] = vNewMetadata
302 vsRowNames = c(vsRowNames,paste(names(dfInput)[i],sLevel,sep="_"))
303 }
304 }
305 }
306
307 # Convert to data frame
308 dfDummyMetadata = as.data.frame(sapply(listMetadata,rbind))
309 names(dfDummyMetadata) = vsRowNames
310 iNumberMetadata = ncol(dfDummyMetadata)
311
312 # Data to use in ordination in NMDS
313 # All cleaned bug data
314 dfData = dfInput[viData]
315
316 # If rotating the ordination by a metadata
317 # 1. Add in the metadata as a bug
318 # 2. Multiply the bug by the weight
319 # 3. Push this through the NMDS
320 if(!is.null(sRotateByMetadata))
321 {
322 vsRotateMetadata = unlist(strsplit(sRotateByMetadata,","))
323 sMetadata = vsRotateMetadata[1]
324 dWeight = as.numeric(vsRotateMetadata[2])
325 sOrdinationMetadata = dfDummyMetadata[sMetadata]*dWeight
326 dfData[sMetadata] = sOrdinationMetadata
327 }
328
329 # Run NMDS on bug data (Default B-C)
330 # Will have species and points because working off of raw data
331 mNMDSData = metaMDS(dfData,k=2)
332
333 ## Make shapes
334 # Defines the shapes and the metadata they are based on
335 # Metadata to use as shapes
336 lShapeInfo = funcMakeShapes(dfInput=dfInput, sShapeBy=sShapeBy, sShapes=sShapes, cDefaultShape=sDefaultMarker)
337
338 sMetadataShape = lShapeInfo[["ID"]]
339 vsShapeValues = lShapeInfo[["Values"]]
340 vsShapeShapes = lShapeInfo[["Shapes"]]
341 vsShapes = lShapeInfo[["PlotShapes"]]
342 cDefaultShape = lShapeInfo[["DefaultShape"]]
343
344 # Colors
345 vsColors = rep(cDefaultColor,nrow(dfInput))
346 vsColorValues = c()
347 vsColorRBG = c()
348 if(!is.null(sColorBy))
349 {
350 vsColorValues = paste(sort(unique(unlist(dfInput[[sColorBy]])),na.last=TRUE))
351 iLengthColorValues = length(vsColorValues)
352
353 vsColorRBG = lapply(1:iLengthColorValues/iLengthColorValues,colorRamp(vsColorRange))
354 vsColorRBG = unlist(lapply(vsColorRBG, function(x) rgb(x[1]/255,x[2]/255,x[3]/255)))
355
356 for(iColor in 1:length(vsColorRBG))
357 {
358 vsColors[which(paste(dfInput[[sColorBy]])==vsColorValues[iColor])]=vsColorRBG[iColor]
359 }
360
361 #If NAs are seperately given color, then color here
362 if(!is.null(sPlotNAColor))
363 {
364 vsColors[which(is.na(dfInput[[sColorBy]]))] = sPlotNAColor
365 vsColorRBG[which(vsColorValues=="NA")] = sPlotNAColor
366 }
367 }
368
369 print("names(dfDummyMetadata)")
370 print(names(dfDummyMetadata))
371
372 # Reduce the bugs down to the ones in the list to be plotted
373 viBugsToPlot = which(row.names(mNMDSData$species) %in% vsBugsToPlot)
374 viMetadataDummy = which(names(dfDummyMetadata) %in% vsMetadata)
375
376 print("viBugsToPlot")
377 print(viBugsToPlot)
378 print("viMetadataDummy")
379 print(names(dfDummyMetadata)[viMetadataDummy])
380
381 # Build the matrix of metadata coordinates
382 mMetadataCoordinates = matrix(rep(NA, iNumberMetadata*2),nrow=iNumberMetadata)
383 for( i in 1:iNumberMetadata )
384 {
385 lxReturn = NA
386 if( i %in% viContinuousMetadata )
387 {
388 lxReturn = funcGetMaximumForMetadatum(dfDummyMetadata[[i]],mNMDSData$points)
389 } else {
390 lxReturn = funcGetCentroidForMetadatum(dfDummyMetadata[[i]],mNMDSData$points)
391 }
392 mMetadataCoordinates[i,] = c(lxReturn$x,lxReturn$y)
393 }
394 row.names(mMetadataCoordinates) = vsRowNames
395
396 # Plot the biplot with the centroid constructed metadata coordinates
397 if(length(viMetadataDummy)==0)
398 {
399 viMetadataDummy = 1:nrow(mMetadataCoordinates)
400 }
401
402 # Plot samples
403 # Make output name
404 if(is.null(sOutputFileName))
405 {
406 viPeriods = which(sInputFileName==".")
407 if(length(viPeriods)>0)
408 {
409 sOutputFileName = paste(OutputFileName[1:viPeriods[length(viPeriods)]],"pdf",sep=".")
410 } else {
411 sOutputFileName = paste(sInputFileName,"pdf",sep=".")
412 }
413 }
414
415 pdf(sOutputFileName, useDingbats=FALSE)
416 plot(mNMDSData$points, xlab=paste("NMDS1","Stress=",mNMDSData$stress), ylab="NMDS2", pch=vsShapes, col=vsColors)
417 title(sTitle,line=3)
418
419 # Plot Bugs
420 mPlotBugs = mNMDSData$species[viBugsToPlot,]
421 if(length(viBugsToPlot)==1)
422 {
423 text(x=mPlotBugs[1],y=mPlotBugs[2],labels=row.names(mNMDSData$species)[viBugsToPlot],col=sTextColor)
424 } else if(length(viBugsToPlot)>1){
425 text(x=mPlotBugs[,1],y=mPlotBugs[,2],labels=row.names(mNMDSData$species)[viBugsToPlot],col=sTextColor)
426 }
427
428 # Add alternative axes
429 axis(3, col=sArrowColor)
430 axis(4, col=sArrowColor)
431 box(col = "black")
432
433 # Plot Metadata
434 if(length(viMetadataDummy)>0)
435 {
436 for(i in viMetadataDummy)
437 {
438 curCoordinates = mMetadataCoordinates[i,]
439 curCoordinates = curCoordinates * dResizeArrow
440 # Plot Arrow
441 arrows(0,0, curCoordinates[1] * 0.8, curCoordinates[2] * 0.8, col=sArrowColor, length=0.1 )
442 }
443 # Plot text
444 if(length(viMetadataDummy)==1)
445 {
446 text(x=mMetadataCoordinates[viMetadataDummy,][1]*dResizeArrow*0.8, y=mMetadataCoordinates[viMetadataDummy,][2]*dResizeArrow*0.8, labels=row.names(mMetadataCoordinates)[viMetadataDummy],col=sArrowTextColor)
447 } else {
448 text(x=mMetadataCoordinates[viMetadataDummy,1]*dResizeArrow*0.8, y=mMetadataCoordinates[viMetadataDummy,2]*dResizeArrow*0.8, labels=row.names(mMetadataCoordinates)[viMetadataDummy],col=sArrowTextColor)
449 }
450 }
451
452 # Create Legend
453 # The text default is the colorMetadata_level (one per level) plus the ShapeMetadata_level (one per level)
454 # The color default is already determined colors plus grey for shapes.
455 sLegendText = c(paste(vsColorValues,sColorBy,sep="_"),paste(sMetadataShape,vsShapeValues,sep="_"))
456 sLegendColors = c(vsColorRBG,rep(cDefaultColor,length(vsShapeValues)))
457
458 # If the color values are numeric
459 # Too many values may be given in the legend (given they may be a continuous range of values)
460 # To reduce this they are summarized instead, given the colors and values for the extreme ends.
461 if( !sum( is.na( as.numeric( vsColorValues[ which( !is.na( vsColorValues ) ) ] ) ) ) )
462 {
463 vdNumericColors = as.numeric( vsColorValues )
464 vdNumericColors = vdNumericColors[ which( !is.na( vdNumericColors ) ) ]
465 vdSortedNumericColors = sort( vdNumericColors )
466 sLegendText = c( paste( sColorBy, vdSortedNumericColors[ 1 ], sep="_" ),
467 paste( sColorBy, vdSortedNumericColors[ length(vdSortedNumericColors) ], sep="_" ),
468 paste( sMetadataShape, vsShapeValues, sep="_" ) )
469 sLegendColors = c(vsColorRBG[ which( vdNumericColors == vdSortedNumericColors[ 1 ] )[ 1 ] ],
470 vsColorRBG[ which( vdNumericColors == vdSortedNumericColors[ length( vdSortedNumericColors ) ] )[ 1 ] ],
471 rep(cDefaultColor,length(vsShapeValues)))
472 }
473 sLegendShapes = c( rep( cDefaultShape, length( sLegendText ) - length( vsShapeShapes ) ), vsShapeShapes )
474
475 # If any legend text was constructed then make the legend.
476 if( length( sLegendText ) >0 )
477 {
478 legend( "topright", legend = sLegendText, pch = sLegendShapes, col = sLegendColors )
479 }
480
481 # Original biplot call if you want to check the custom plotting of the script
482 # There will be one difference where the biplot call scales an axis, this one does not. In relation to the axes, the points, text and arrows should still match.
483 # Axes to the top and right are for the arrow, others are for markers and bug names.
484 #biplot(mNMDSData$points,mMetadataCoordinates[viMetadataDummy,],xlabs=vsShapes,xlab=paste("MDS1","Stress=",mNMDSData$stress),main="Biplot function Bugs and Sampes - Metadata Plotted with Centroids")
485 dev.off()
486 }
487
488 # This is the equivalent of __name__ == "__main__" in Python.
489 # That is, if it's true we're being called as a command line script;
490 # if it's false, we're being sourced or otherwise included, such as for
491 # library or inlinedocs.
492 if( identical( environment( ), globalenv( ) ) &&
493 !length( grep( "^source\\(", sys.calls( ) ) ) )
494 {
495 lsArgs <- parse_args( pArgs, positional_arguments=TRUE )
496
497 funcDoBiplot(
498 sBugs = lsArgs$options$sBugs,
499 sMetadata = lsArgs$options$sMetadata,
500 sColorBy = lsArgs$options$sColorBy,
501 sColorRange = lsArgs$options$sColorRange,
502 sTextColor = lsArgs$options$sTextColor,
503 sArrowColor = lsArgs$options$sArrowColor,
504 sArrowTextColor = lsArgs$options$sArrowTextColor,
505 sPlotNAColor = lsArgs$options$sPlotNAColor,
506 sShapeBy = lsArgs$options$sShapeBy,
507 sShapes = lsArgs$options$sShapes,
508 sDefaultMarker = lsArgs$options$sDefaultMarker,
509 sRotateByMetadata = lsArgs$options$sRotateByMetadata,
510 dResizeArrow = lsArgs$options$dResizeArrow,
511 sTitle = lsArgs$options$sTitle,
512 sInputFileName = lsArgs$args[2],
513 sLastMetadata = lsArgs$args[1],
514 sOutputFileName = lsArgs$options$sOutputFileName)
515 }