## ----options, include=FALSE, echo=FALSE------------------------------------
library(BiocStyle)
knitr::opts_chunk$set(warning=FALSE, error=FALSE, message=FALSE)

## ----construct-------------------------------------------------------------
library(SingleCellExperiment)
counts <- matrix(rpois(100, lambda = 10), ncol=10, nrow=10)
sce <- SingleCellExperiment(assays = list(counts = counts))
sce

## ----coerce----------------------------------------------------------------
se <- SummarizedExperiment(list(counts=counts))
as(se, "SingleCellExperiment")

## ----fluidigm--------------------------------------------------------------
library(scRNAseq)
sce <- ReprocessedAllenData("tophat_counts")
sce

## ----subset----------------------------------------------------------------
counts <- assay(sce, "tophat_counts")
libsizes <- colSums(counts)
size.factors <- libsizes/mean(libsizes)
logcounts(sce) <- log2(t(t(counts)/size.factors) + 1)
assayNames(sce)

## ----pca-------------------------------------------------------------------
pca_data <- prcomp(t(logcounts(sce)), rank=50)

library(Rtsne)
set.seed(5252)
tsne_data <- Rtsne(pca_data$x[,1:50], pca = FALSE)

reducedDims(sce) <- list(PCA=pca_data$x, TSNE=tsne_data$Y)
sce

## --------------------------------------------------------------------------
reducedDims(sce)
reducedDimNames(sce)
head(reducedDim(sce, "PCA")[,1:2])
head(reducedDim(sce, "TSNE")[,1:2])

## --------------------------------------------------------------------------
dim(reducedDim(sce, "PCA"))
dim(reducedDim(sce[,1:10], "PCA"))

## --------------------------------------------------------------------------
counts(sce) <- assay(sce, "tophat_counts")
sce
dim(counts(sce))

## --------------------------------------------------------------------------
is.spike <- grepl("^ERCC-", rownames(sce))
sce <- splitAltExps(sce, ifelse(is.spike, "ERCC", "gene"))
altExpNames(sce)

## --------------------------------------------------------------------------
altExp(sce)

## --------------------------------------------------------------------------
rowData(altExp(sce))$concentration <- runif(nrow(altExp(sce)))
rowData(altExp(sce))
rowData(sce)

## --------------------------------------------------------------------------
sessionInfo()