read.FCS {flowCore}R Documentation

Read an FCS file

Description

Check validity and Read Data File Standard for Flow Cytometry

Usage


   isFCSfile(files)

   read.FCS(filename, transformation="linearize", which.lines=NULL,
            debug=FALSE, alter.names=FALSE, column.pattern=NULL,
            decades=0, ncdf = FALSE, min.limit=NULL)

   cleanup()

Arguments

files A vector of filenames
filename Character of length 1: filename
transformation An character string that defines the type of transformation. Valid values are linearize (default) or scale.The linearize tramsformation applies the appropriate power transform to the data while the scale transformation scales all columns to $[0,10^decades]$. defaulting to decades=0 as in the FCS4 specification. A logical can also be used: TRUE is equal to linearize and FALSE correponds to no transformation.
which.lines Numeric vector to specify the indices of the lines to be read. If NULL all the records are read, if of length 1, a random sample of the size indicated by which.lines is read in.
debug boolean indicating whether or not to print the debugging statements, default is TRUE
alter.names boolean indicating whether or not we should rename the columns to valid R names using make.names. The default is FALSE.
column.pattern An optional regular expression defining parameters we should keep when loading the file. The default is NULL.
decades When scaling is activated, the number of decades to use for the output.
ncdf Instead of reading all data into memory, this switches to file-based data storage. A netCDF file is creates for each flowFrame in the .flowCoreNcdf subdirectory. For large data sets this significantly reduces the memory profile of the R session, to the cost of speed and disk space. The exprs and exprs<- methods make sure that the user always gets a matrix of data values. Please note that currently all operations that call exprs<-, either explicitely or implicitely, will result in the creation of a new netCDF file. This behaviour may change in the future. Currently the software does not remove any of the netCDF files and it is up to the user to do clean up. The easiest way to do that is to delete the whole netCDF directory. To this end, one can envoke the cleanup function.
min.limit The minimum value in the data range that is allowed. Some instruments produce extreme artifactual values. The positive data range for each parameter is completely defined by the measurement range of the instrument and all larger values are set to this threshold. The lower data boundary is not that well defined, since compensation might shift some values below the original measurement range of the instrument. The default value of -111 copies the behavior of flowJo. It can be set to an arbitray number or to NULL, in which case the original values are kept.

Details

The function isFCSfile determines whether its arguments are valid FCS files.

The function read.FCS works with the output of the FACS machine software from a number of vendors (FCS 2.0, FCS 3.0 and List Mode Data LMD). However, the FCS 3.0 standard includes some options that are not yet implemented in this function. If you need extensions, please let me know. The output of the function is an object of class flowFrame.

For specifications of FCS 3.0 see http://www.isac-net.org and the file ../doc/fcs3.html in the doc directory of the package.

The nlines and sampling arguments allow you to read a subset of the record as you might not want to read the thousands of events recorded in the FCS file.

The which.lines argument allows you to read a specific number of records.

Value

isFCSfile returns a logical vector.
read.FCS returns an object of class flowFrame that contains the data in the exprs slot, the parameters monitored in the parameters slot and the keywords and value saved in the header of the FCS file.

Author(s)

F. Hahne, N.Le Meur

See Also

link[flowCore]{read.flowSet}

Examples

## a sample file
fcsFile <- system.file("extdata", "0877408774.B08", package="flowCore")

## read file and linearize values
samp <-  read.FCS(fcsFile, transformation="linearize")
exprs(samp[1:3,])
description(samp)[3:6]
class(samp)

## Only read in lines 2 to 5
subset <- read.FCS(fcsFile, which.lines=2:5, transformation="linearize")
exprs(subset)

## Read in a random sample of 100 lines
subset <- read.FCS(fcsFile, which.lines=100, transformation="linearize")
nrow(subset)


[Package flowCore version 1.10.0 Index]