# 6. Gene Coexpression

Code for how to extract the cells that co-express two genes.

1. Get the ensembl id for your genes
``````rowData(sce)[which(rowData(sce)[,2] %in% "Lefty1"),1]
rowData(sce)[which(rowData(sce)[,2] %in% "Lefty2"),1]
``````
1. Which rows of the sce database are those genes found in?
``````which(rowData(sce)[,2] %in% "Lefty1")
which(rowData(sce)[,2] %in% "Lefty2")
# double check you are right
rowData(sce)[1571:1573,]
``````
1. Create a new sce object with just those 2 genes.
column 1 is Lefty1 and column 2 is Lefty2
``````leftys <- sce[c(1573,1571),]
``````
1. Binary Approach: Double vs not double
``````keepers <- counts(leftys)[1,] != 0 & counts(leftys)[2,] != 0
table(keepers)
# TRUE = double, FALSE = not double

# create an empty vector with enough spots for all the cells
lefty.colors <- c(1:length(keepers))

# run a for loop, this is an advanced topic,
# purpose: Convert TRUE to "double" and FALSE to "not double".
for (i in c(1:length(keepers))){
if (keepers[i]==TRUE){lefty.colors[i]='double'}
else {lefty.colors[i]='not double'}
}

# Add your new cell label back to the sce object to plot altogether
colData(sce)\$Leftys <- lefty.colors
plotReducedDim(object = sce, dimred="umap", colour_by="Leftys")
``````
1. 4-way approach: Double, Singles, Zeros
``````# create an empty vector with enough spots for all the cells
coexpress <- as.data.frame(c(1:ncol(sce)))

# run a for loop, this is an advanced topic,
# purpose: Determine if cells double positive, single positive, or neither.
for (i in c(1:ncol(sce))){
if (counts(leftys)[1,i] != 0 & counts(leftys)[2,i] != 0){
coexpress\$Group[i]="double" }
else if (counts(leftys)[1,i] > 0 & counts(leftys)[2,i] == 0){
coexpress\$Group[i]="Lefty1"}
else if (counts(leftys)[1,i] == 0 & counts(leftys)[2,i] > 0){
coexpress\$Group[i]="Lefty2"}
else {
coexpress\$Group[i]="neither"}
}

# check how many cells fall in each category
table(coexpress\$Group)

# Add it to the sce cell data then plot
colData(sce)\$LeftyGroup <- coexpress\$Group
plotReducedDim(object = sce, dimred="umap", colour_by="Group")
``````

Note: 4-way may need some editing.

1 Like