This combination of parameters throws an error. We need to determine what should be returned (if anything) and fix it, or report a better message.

The "continuous" portion of the object isn't being initialized correctly.

Any idea what could cause this error?

assay.f <- SingleCellAssay::filter(assay, groups="time")
Error in apply(abs(z.exprs) > sigmaContinuous, 1, sum, na.rm = T) :
dim(X) must have a positive length

the assay was built like:

assay <- SingleCellAssay(dataframe=ncounts.m, idvars="cell", cellvars="time", primerid="tx", featurevars="symbol", measurement="log.count")

time is a factor with 5 levels of which 3 are used.

Thanks

We'd like to be able to do things like:

O is a SingleCellAssay object.

mapping(O)[["primerid"]]<-"foo"

which would retrieve the mapping from O, access the "primerid" mapping and set it to "foo", then update the @mapping slot in O.
Needs a replace method for mapping()<-
and [[<-

No need to change the mapping class, just write a couple of extra methods.

Can you check out why the devel branch is buggy? Given that the variance shrinkage is here, it will eventually be merged into master. I can't sort out what problem is causing roxygen to fail, ie. whether it's the code or a bug in roxygen. But there are a lot of roxyen like comments peppering the code which may be causing problems.

Hi

I've just started using SingleCellAssay, and have found that loading the library after loading reshape2 broke my calls to melt() in functions that I had written (session info below).

I think it's a bit impolite to mask reshape2's melt(), which is a much more generic and widely-used function. Perhaps you could tweak melt.SingleCellAssay so it plays more nicely with reshape2?

Best wishes
Davis

sessionInfo()
R version 3.1.1 (2014-07-10)
Platform: x86_64-apple-darwin13.1.0 (64-bit)

locale:
[1] C

attached base packages:
[1] parallel grid stats graphics grDevices utils datasets
[8] methods base

other attached packages:
[1] SingleCellAssay_0.861 scde_1.2.0 flexmix_2.3-12
[4] lattice_0.20-29 gridExtra_0.9.1 RankProd_2.38.0
[7] devtools_1.6.1 ggthemes_1.7.0 beeswarm_0.1.6
[10] magrittr_1.0.1 scater_0.1 roxygen2_4.0.2
[13] knitr_1.7.8 wesanderson_0.3 RColorBrewer_1.0-5
[16] edgeR_3.8.2 limma_3.22.1 dplyr_0.3.0.2
[19] lubridate_1.3.3 stringr_0.6.2 ggplot2_1.0.0
[22] plyr_1.8.1 reshape2_1.4 ProjectTemplate_0.6

loaded via a namespace (and not attached):
[1] Biobase_2.26.0 BiocGenerics_0.12.0 Cairo_1.5-6
[4] DBI_0.3.1 MASS_7.3-35 Rcpp_0.11.3
[7] RcppArmadillo_0.4.500.0 Rook_1.1-1 SparseM_1.05
[10] abind_1.4-0 assertthat_0.1 brew_1.0-6
[13] car_2.0-21 chron_2.3-45 colorspace_1.2-4
[16] compiler_3.1.1 data.table_1.9.4 digest_0.6.4
[19] evaluate_0.5.5 formatR_1.0 gdata_2.13.3
[22] gtable_0.1.2 gtools_3.4.1 labeling_0.3
[25] latticeExtra_0.6-26 lazyeval_0.1.9 memoise_0.2.1
[28] modeltools_0.2-21 munsell_0.4.2 nnet_7.3-8
[31] proto_0.3-10 quantreg_5.05 reshape_0.8.5
[34] rjson_0.2.15 scales_0.2.4 stats4_3.1.1
[37] tools_3.1.1

It seems the SCA constructor depends on the order of rows in the input data.frame since the internal representation has been changed to a data.table

Test script:

data(vbeta)
vbeta <- computeEtFromCt(vbeta)
# okay when data is ordered by wells, primers
vbeta.sca <- SingleCellAssay(vbeta, idvars = c("Subject.ID", "Chip.Number", "Well"),
primerid = "Gene", measurement = "Et", geneid = "Gene",
cellvars = c("Number.of.Cells", "Population"), phenovars = c("Stim.Condition",
"Time"), id = "vbeta all")
print(vbeta.sca)
# fails after data is ordered by primers,wells
reordering<-order(vbeta[,c("Gene","Subject.ID", "Chip.Number", "Well")])
vbeta.reordered=vbeta[reordering,]
vbeta.reord.sca <- SingleCellAssay(vbeta.reordered, idvars = c("Subject.ID", "Chip.Number", "Well"),
primerid = "Gene", measurement = "Et", geneid = "Gene",
cellvars = c("Number.of.Cells", "Population"), phenovars = c("Stim.Condition",
"Time"), id = "vbeta all reordered")

Output:

> data(vbeta)
> vbeta <- computeEtFromCt(vbeta)
> # okay when data is ordered by wells, primers
> vbeta.sca <- SingleCellAssay(vbeta, idvars = c("Subject.ID", "Chip.Number", "Well"),
+ primerid = "Gene", measurement = "Et", geneid = "Gene",
+ cellvars = c("Number.of.Cells", "Population"), phenovars = c("Stim.Condition",
+ "Time"), id = "vbeta all")
> print(vbeta.sca)
SingleCellAssay  id:  vbeta all 
 456  wells;  75  features

> # fails after data is ordered by primers,wells
> reordering<-order(vbeta[,c("Gene","Subject.ID", "Chip.Number", "Well")])
> vbeta.reordered=vbeta[reordering,]
> vbeta.reord.sca <- SingleCellAssay(vbeta.reordered, idvars = c("Subject.ID", "Chip.Number", "Well"),
+ primerid = "Gene", measurement = "Et", geneid = "Gene",
+ cellvars = c("Number.of.Cells", "Population"), phenovars = c("Stim.Condition",
+ "Time"), id = "vbeta all")
Error in SingleCellAssay(vbeta.reordered, idvars = c("Subject.ID", "Chip.Number",  : 
  'featurevars' must be keyed by 'primerid'

*I think the function should work on log transformed data. It takes and returns data on the natural scale, but many of the arguments (eg qt or absolute_min) are in terms of the log scale. So we end up just logging the data twice.

*Or we could make log-transforming that data an optional argument if we think that people would ever not want to log-transform

*There's a bug somewhere in lines 166-176 that is causing it not to remove an empty bin if it's the first bin in the sequence. When this happens, the call to density barfs.

*It needs test coverage.

Currently we can import from a matrix without reshaping by starting from an ExpressionSet, but this seems to leave some important attributes unset (like primerid).

Replace method for cData()<- is needed to replace the phenoData assiged to the cData slot

I noticed this package is only available as a .zip file. Because of certain restrictions at work, I am unable to download the package via R and must download and install manually. We use macs here at work and was wondering if this package is available as .tgz or .tar.gz file since every time I try to use R CMD INSTALL in terminal it gives me "non-zero exit status" error. Thanks

I think it would be worthwhile supporting limma-like interface where the model is fit once, and contarsts can be extracted via function like contrasts.fit in limma.
@amcdavid I'm filing this as an enhancement, we can discuss whether and how this should be done.

Currently SingleCellAssay class contains slots of a cmap and fmap, but these are only used under construction when we pass in a data.frame. These probably don't need to be part of the API, but we should think about this a little bit.

thresholdNanoString currently can output a list of diagnostics with the clustering. This should be wrapped into a class.

I've got plotting code that should be turned into a method that dispatches on the thresholdNanoString class.

Issue 1

When calling getConcordance(FluidigmAssay,'ncell'):
Error in UseMethod("melt", data) :
no applicable method for 'melt' applied to an object of class "c('cast_df', 'data.frame')"

from line 115 in Fluidigm-methods.R:

castL[[i]] <- cast(melt(tmp), secondForm, fun.aggregate=fun.cycle)

If I'm not completely wrong this line should be corrected to:

castL[[i]] <- cast(tmp, secondForm, fun.aggregate=fun.cycle)

since tmp is a data.frame and not an SingleCellAssay/FluidigmAssay object.

Issue 2

In case of not replacing NAs by 0s in the expression measurements [exprs(sca)], getConcordance on a FluidigmAssay will return n.exp=0, et=NA for all (gene x experimental unit) combinations that have even a single non-expressed well (and thus in this case a single NA value).
The above problem could be fixed by letting the function expavg from Fluidigm-methods.R remove NAs when summing up. e.g. like this:

original (line 14 in Fluidigm-methods.R)

expavg <- function(x) mean(2^x-1)

fixed (could certainly be done more elegantly):

expavg <- function(x) {res<-mean(2^x-1,na.rm=T);ifelse(!is.finite(res),NA,res)}

I've got some single-cell RNA-seq counts (normalised with DESeq2) which I'd like to model with your mixture of continous and discrete components. I'm having some problems getting any sensible p-values from zlm.SingleCellAssay and I'm not finding the intro vignette or the package documentation particularly helpful. Is there any chance you could put together a small example, as you did for the fluidigm assay?

Thanks,
John.

This works:
nsa <- NanoStringAssay(files.key, paste(pathToExp, 'h9key.csv', sep=''), idvars=c('ID', 'CartridgeID'), primerid='Name', measurement='Count', geneid='Name', featurevars='GeneID', ncells='ncells', cellvars=c('BindingDensity', 'Date', 'FovCount', 'FovCounted', 'plate', 'cycle', 'cellline'), post.process.function=ppf)

This causes many bad things:
nsa <- NanoStringAssay(files.key, paste(pathToExp, 'h9key.csv', sep=''), idvars=c('ID', 'CartridgeID'), primerid='GeneID', measurement='Count', geneid='Name', featurevars='GeneID', ncells='ncells', cellvars=c('BindingDensity', 'Date', 'FovCount', 'FovCounted', 'plate', 'cycle', 'cellline'), post.process.function=ppf)

this induces at least three bugs, none of which I have time to isolate at the moment.

1. mkunique line 407 should be defined unconditionally, not in the if-clause where its presently located
2. But after fixing this, it's apparent that line 245 in AllClasses.R should read if(!all(object@mapNames%in%getMapNames(getMapping(object))))
3. But even after making the fix in 2, things still explode due to validObject not holding.

I am not going to try to debug this further without a small reproducible example that fails.

Implement reader functions for NanoString lanes (RCC files).

Hi guys, I enjoyed seeing you at CYTO. I tried installing SingleCellAssay using install_github('SingleCellAssay', 'RGLab') on Rstudio 0.97.336 and R version 2.15.3. I got an error relating to the vignette compilation. I installed roxygen2, but that didn't help.

I suspect Sweave's syntax has changed, and your "echo=" arguments are no longer allowed. I was able to resolve the error and install from source through the following steps:

1. Download the .zip
2. Open SingleCellAssay_master/vignettes/SingleCellAssay-intro.Rnw
3. Set all "echo=" to FALSE. This affects lines 66, 142, and 213.
4. Install knitr and plyr packages from CRAN just in case
5. Install SingleCellAssay_master package from source

That may not be the best solution, but hopefully it will point you in the right direction to make installation user-friendly again.

Cheers,

Erin

PS: Console output of the error is below.

> install_github('SingleCellAssay', 'RGLab')
Installing github repo(s) SingleCellAssay/master from RGLab
Installing SingleCellAssay.zip from https://github.com/RGLab/SingleCellAssay/archive/master.zip
Installing SingleCellAssay
Installing dependencies for SingleCellAssay:
data.table
Installing package(s) into ‘/Users/esimonds/Library/R/2.15/library’
(as ‘lib’ is unspecified)

  There is a binary version available (and will be installed) but the source version is later:
           binary source
data.table  1.8.8  1.9.2

trying URL 'http://cran.rstudio.com/bin/macosx/leopard/contrib/2.15/data.table_1.8.8.tgz'
Content type 'application/x-gzip' length 1188814 bytes (1.1 Mb)
opened URL
==================================================
downloaded 1.1 Mb


The downloaded binary packages are in
    /var/folders/cm/8hg79fh52fj5w_h_mmk580tc0000gn/T//Rtmpvm5P3P/downloaded_packages
'/Library/Frameworks/R.framework/Resources/bin/R' --vanilla CMD build  \
  '/private/var/folders/cm/8hg79fh52fj5w_h_mmk580tc0000gn/T/Rtmpvm5P3P/SingleCellAssay-master'  \
  --no-manual --no-resave-data 

* checking for file '/private/var/folders/cm/8hg79fh52fj5w_h_mmk580tc0000gn/T/Rtmpvm5P3P/SingleCellAssay-master/DESCRIPTION' ... OK
* preparing 'SingleCellAssay':
* checking DESCRIPTION meta-information ... OK
* installing the package to re-build vignettes
* creating vignettes ... ERROR
Error: processing vignette 'SingleCellAssay-intro.Rnw' failed with diagnostics:
parse error or empty option in
long-example,warning=FALSE, echo=-c(1,2,3)
Execution halted
Error: Command failed (1)

When coefficient is not estimated NA is produced as the result of logFC or getLogFC.

I get errors when trying to run likelihood ratio tests on FluidigmAssays in the current version of the package.

data(vbeta)
vbeta <- computeEtFromCt(vbeta)
vbeta.fa <- FluidigmAssay(vbeta, idvars = c("Subject.ID", "Chip.Number", "Well"),
  primerid = "Gene", measurement = "Et", ncells = "Number.of.Cells", geneid = "Gene",
  cellvars = c("Number.of.Cells", "Population"), phenovars = c("Stim.Condition",
  "Time"), id = "vbeta all")
head(cellData(vbeta.fa)@data)
#            Number.of.Cells            Population Subject.ID Chip.Number Well
#Sub01:1:A01               1 CD154+VbetaResponsive      Sub01           1  A01
#Sub01:1:A02               1 CD154+VbetaResponsive      Sub01           1  A02
#Sub01:1:A03               1 CD154+VbetaResponsive      Sub01           1  A03
#Sub01:1:A04               1 CD154+VbetaResponsive      Sub01           1  A04
#Sub01:1:A05               1 CD154+VbetaResponsive      Sub01           1  A05
#Sub01:1:A06               1 CD154+VbetaResponsive      Sub01           1  A06
#            Stim.Condition Time
#Sub01:1:A01      Stim(SEB)   12
#Sub01:1:A02      Stim(SEB)   12
#Sub01:1:A03      Stim(SEB)   12
#Sub01:1:A04      Stim(SEB)   12
#Sub01:1:A05      Stim(SEB)   12
#Sub01:1:A06      Stim(SEB)   12
LRT(vbeta.fa,comp='Stim.Condition',ref='Stim(SEB)')
#Error in split(x[[measure]], x$pheno.order, drop = FALSE) [line 130 in lrtest.R]: 
#  error in evaluating the argument 'x' in selecting a method for function 'split': Error in x[[measure]] : subscript out of bounds

If x is a SingleCellAssay object, then x[,y] coerces y improperly.
Factors are being coerced into integer indexes. This causes confusing behavior. We should throw an error rather than coerce.

Also, check to see why fData(x)$primerid is getting promoted to a factor in the first place.

Will profile it to try to find the bottleneck.

R version 3.2
SingleCellAssay 0.78
Bioconductor 2.14 (development)
combine is failing within filter.

A lot of reproducible code is no longer reproducible. I'm now unable to deliver an analysis on time.

>lapp

$`High-DCVax-001:154- CM`
FluidigmAssay  on layer  Et 
 1  Layers;  16  wells;  96  features
 id:  High-DCVax-001:154- CM 

$`High-DCVax-001:154+ CM`
FluidigmAssay  on layer  Et 
 1  Layers;  256  wells;  96  features
 id:  High-DCVax-001:154+ CM 

$`Low-DCVax-001:154+ CM`
FluidigmAssay  on layer  Et 
 1  Layers;  100  wells;  96  features
 id:  Low-DCVax-001:154+ CM 

$`Mid-DCVax-001:154+ CM`
FluidigmAssay  on layer  Et 
 1  Layers;  119  wells;  96  features
 id:  Mid-DCVax-001:154+ CM 

$`Placebo:154+ CM`
FluidigmAssay  on layer  Et 
 1  Layers;  90  wells;  96  features
 id:  Placebo:154+ CM 

combine(as(lapp, "SCASet"))

Error in combine(as(lapp, "SCASet")) : 
  error in evaluating the argument 'x' in selecting a method for function 'combine': Error in asMethod(object) : 
  All members of 'x' must inherit from 'SingleCellAssay'

Document arguments with ?? in 790f88f

The new data.table constructor has a bunch of issues in the tests (the constructed object is essentially gibberish right now) so I moved that code to devel and reverted master.

Agreed that we need to make RNAseq construction faster. Would writing a constructor that takes an expression matrix and fData and cData be a good solution?

Analysis functions should either 1) take na.rm, or 2) we should have a slot in SingleCellAssay objects that specify how NAs will be treated in analysis.

I'd lean towards 1, because adding a slot is likely to be a bit invasive at this point. @gfinak your thoughts?

And either use glm.fit2 to ensure convergence or test for decrease in deviance after fitting

Need to write S3 generic for melt if one does not already exist when S4 generic is defined.

doubleid <- data.frame(id1=1:3, id2=1:3, et=rep(3, 3), f1=rep('A', 3))
smallsc <- SingleCellAssay(doubleid, idvars='id1', geneid='f1', primerid='f1', measurement='et')
spl <- split(smallsc, 'id1')
c2 <- combine(spl[[1]], spl[[2]], spl[[3]]) #Works fine

c2 <- do.call(combine, spl@set)

cannot coerce type 'closure' to vector of type 'character'
1: do.call(combine, spl@set)
2: structure(function (x, y, ...)
{
if (length(list(...)) > 0L) {
callGeneric(x, do.call(callGeneric, list(y, ...)))
}
else {
standardGeneric("combine")
}
}, generic = structure("combine", package = "BiocGenerics"), package = "BiocGenerics", group = list(), valueClass = character(0), signature = c("x", "y"), default = \001NULL\001, skeleton = function (x, y, ...)
stop("invalid call in method dispatch to "combine" (no default method)", domain = NA)(x, y, ...), class = structure("nonstandardGenericFunction", package = "methods"))(1 = <S4 object of class structure("SingleCellAssay", package = "SingleCellAssay")>, 2 = <S4 object of class structure("SingleCellAssay", package = "SingleCellAssay")>, 3 = <S4 object of class structure("SingleCellAssay", package = "SingleCellAssay")>)
3: callGeneric(x, do.call(callGeneric, list(y, ...)))
4: as.character(call[[1L]])
5: as.character.default(call[[1L]])

1. More complicated zlm examples
2. RNAseq examples

Returns data out of order.

The logical matrix returned when apply_filter=FALSE and filt_control(filter=TRUE) is labelled by wellKey. This is not very legible. We should label it with the set of idvars before feeding it to the user.

The package currently fails regression tests.

During R CMD build, I get

Error in setMethod("fit", signature = c(object = "GLMlike", response = "missing"),  :
  no existing definition for function ‘fit’
Error : unable to load R code in package ‘SingleCellAssay’
ERROR: lazy loading failed for package ‘SingleCellAssay’

Create method for NanostringAssay objects that updates the "measurement" field of the mapping to contain thresholded data.

Signature:
threshold(NanoObject, groups, manualThresholds=NULL)
NanoObject: NanostringAssay object
groups: optional grouping variable, thresholds will be calculating per each group
manualThresholds: dataframe giving manual cutoffs for various genes/groups. If not provided gaussian mixture models are used instead.

Slots/Mappings required:
-Number.of.Cells (already available if we subclass FluidigmAssay)
-raw

Questions:
-Do we want to mark indeterminate values as NA? (values without definitive cluster membership posterior probabilities)
-We should log2 +1 transform the data either here or when the NanostringAssay object is constructed

I'll start by constructing a FluidigmAssay from NanoString data, then work on specializing a NanoStringAssay class that captures anything which is missing.

Okay, I'm able to construct a FluidigmAssay object without much trouble.

-Test coverage
-Needs to be profiled
-Support arbitrary contrasts
-Needs friendlier interface
-Look for a way to avoid the bootstrap
-Return average correlation in test set

When using combine to add to the featureData or cellData, we should make it harder to incorrectly add data by matching by row.names or a common field if possible.

New method for fData replacement would be consistent with cData.

constructors were expected to take certain named parameters so that combine would work properly and construct the appropriate type of assay. I broke this with the NanostringAssay constructor (bacause I forgot about it).
Needs to be fixed.

Lucas mentioned that some genes were not unique and should be summed to get a count within the lane.
Calling exprs() on a NanoStringAssay object tries to construct a matrix of columns = number of genes (or maybe rows). This does not agree with the unique number of genes and errors out.

We should either sum over the replicates before constructing the object, or define a unique primer id , or change the way exprs() gets the number of rows to take replicate genes into account.

We probably should submit a patch to car for lht rather than fixing this here.

But on the bright side, now the model fit objects from lmer support drop1. so our LRT testing code should work out of the box.