Figure 5 - Table S5D

Author

Ashok Patowary

isoT4k_filter = readRDS("data/working/isoT4k_filter.rds")
levels(isoT4k_filter$orig.ident) = gsub("\\.", "-", levels(isoT4k_filter$orig.ident)) # ExM.U becomes ExM-U
Loading required package: SeuratObject
Attaching sp
Loading required package: Seurat
isoT4k_filter
An object of class Seurat 
141112 features across 4167 samples within 2 assays 
Active assay: RNA (87162 features, 0 variable features)
 1 other assay present: SCT
 2 dimensional reductions calculated: pca, umap
################Identifying gene with different isoforms significantly expressed in different cell types
#####Code is adapted from https://doi.org/10.1038/s41467-020-17800-6


clusters <- c("End", "ExDp1", "ExDp2", "ExM", "ExM.U", "ExN", "InCGE", "InMGE", "IP", "Mic", "OPC", "oRG", "Per", "PgG2M", "PgS", "vRG")
Idents(isoT4k_filter) <- isoT4k_filter$CellType

total <- data.frame()
totaladj <- data.frame()
for (i in (1:15)){
  k <- i+1
  for (j in (k:16)){
    if(i != j){
      print(paste(i, " ", j, " ",clusters[i], " vs ", clusters[j], sep=""))
      
      markers <- FindMarkers(object = isoT4k_filter , ident.1=clusters[i], ident.2=clusters[j])
      markers$cluster <- clusters[j]
      markers$contrast <- paste(clusters[i], "vs", clusters[j], sep=" ")
      markers[which(markers$avg_log2FC>0),]$cluster <- clusters[i]
      markers$geneId <- sapply(strsplit(rownames(markers), "-ENS|-TAL"), `[`, 1)
      markers$transcriptId <- sapply(strsplit(rownames(markers), "-ENS|-TAL"), `[`, 2)
      markers <- markers[markers$p_val < 0.05,]
      all.genes <- unique(markers$geneId)
      for (k in (1:length(all.genes))){
        sub <- markers[which(markers$geneId == all.genes[k]),]
        nb.clusters <- unique(sub$cluster)
        nb.transcripts <- unique(sub$transcriptId)
        
        if(length(nb.clusters) > 1 & length(nb.transcripts) > 1){
          total <- rbind(total, sub)
        }
      }
      
      markers <- markers[markers$p_val_adj < 0.05,]
      all.genes <- unique(markers$geneId)
      for (k in (1:length(all.genes))){
        sub <- markers[which(markers$geneId == all.genes[k]),]
        nb.clusters <- unique(sub$cluster)
        nb.transcripts <- unique(sub$transcriptId)
        
        if(length(nb.clusters) > 1 & length(nb.transcripts) > 1){
          totaladj <- rbind(totaladj, sub)
        }
      }
      print (dim(total))
      print (length(unique(total$geneId)))
    }
  }
}
[1] "1 2 End vs ExDp1"
[1] 10  9
[1] 4
[1] "1 3 End vs ExDp2"
[1] 10  9
[1] 4
[1] "1 4 End vs ExM"
[1] 25  9
[1] 7
[1] "1 5 End vs ExM.U"
[1] 44  9
[1] 10
[1] "1 6 End vs ExN"
[1] 56  9
[1] 10
[1] "1 7 End vs InCGE"
[1] 66  9
[1] 10
[1] "1 8 End vs InMGE"
[1] 75  9
[1] 10
[1] "1 9 End vs IP"
[1] 92  9
[1] 10
[1] "1 10 End vs Mic"
[1] 92  9
[1] 10
[1] "1 11 End vs OPC"
[1] 94  9
[1] 10
[1] "1 12 End vs oRG"
[1] 110   9
[1] 11
[1] "1 13 End vs Per"
[1] 110   9
[1] 11
[1] "1 14 End vs PgG2M"
[1] 135   9
[1] 13
[1] "1 15 End vs PgS"
[1] 155   9
[1] 14
[1] "1 16 End vs vRG"
[1] 172   9
[1] 15
[1] "2 3 ExDp1 vs ExDp2"
[1] 172   9
[1] 15
[1] "2 4 ExDp1 vs ExM"
[1] 178   9
[1] 18
[1] "2 5 ExDp1 vs ExM.U"
[1] 180   9
[1] 19
[1] "2 6 ExDp1 vs ExN"
[1] 184   9
[1] 20
[1] "2 7 ExDp1 vs InCGE"
[1] 188   9
[1] 22
[1] "2 8 ExDp1 vs InMGE"
[1] 192   9
[1] 24
[1] "2 9 ExDp1 vs IP"
[1] 200   9
[1] 24
[1] "2 10 ExDp1 vs Mic"
[1] 221   9
[1] 29
[1] "2 11 ExDp1 vs OPC"
[1] 221   9
[1] 29
[1] "2 12 ExDp1 vs oRG"
[1] 246   9
[1] 38
[1] "2 13 ExDp1 vs Per"
[1] 256   9
[1] 40
[1] "2 14 ExDp1 vs PgG2M"
[1] 273   9
[1] 41
[1] "2 15 ExDp1 vs PgS"
[1] 283   9
[1] 41
[1] "2 16 ExDp1 vs vRG"
[1] 302   9
[1] 41
[1] "3 4 ExDp2 vs ExM"
[1] 306   9
[1] 43
[1] "3 5 ExDp2 vs ExM.U"
[1] 314   9
[1] 45
[1] "3 6 ExDp2 vs ExN"
[1] 318   9
[1] 45
[1] "3 7 ExDp2 vs InCGE"
[1] 320   9
[1] 46
[1] "3 8 ExDp2 vs InMGE"
[1] 322   9
[1] 47
[1] "3 9 ExDp2 vs IP"
[1] 330   9
[1] 48
[1] "3 10 ExDp2 vs Mic"
[1] 332   9
[1] 48
[1] "3 11 ExDp2 vs OPC"
[1] 334   9
[1] 48
[1] "3 12 ExDp2 vs oRG"
[1] 346   9
[1] 48
[1] "3 13 ExDp2 vs Per"
[1] 346   9
[1] 48
[1] "3 14 ExDp2 vs PgG2M"
[1] 350   9
[1] 49
[1] "3 15 ExDp2 vs PgS"
[1] 363   9
[1] 50
[1] "3 16 ExDp2 vs vRG"
[1] 374   9
[1] 50
[1] "4 5 ExM vs ExM.U"
[1] 376   9
[1] 50
[1] "4 6 ExM vs ExN"
[1] 378   9
[1] 51
[1] "4 7 ExM vs InCGE"
[1] 378   9
[1] 51
[1] "4 8 ExM vs InMGE"
[1] 383   9
[1] 52
[1] "4 9 ExM vs IP"
[1] 387   9
[1] 53
[1] "4 10 ExM vs Mic"
[1] 405   9
[1] 54
[1] "4 11 ExM vs OPC"
[1] 411   9
[1] 56
[1] "4 12 ExM vs oRG"
[1] 434   9
[1] 59
[1] "4 13 ExM vs Per"
[1] 447   9
[1] 60
[1] "4 14 ExM vs PgG2M"
[1] 455   9
[1] 60
[1] "4 15 ExM vs PgS"
[1] 460   9
[1] 60
[1] "4 16 ExM vs vRG"
[1] 478   9
[1] 61
[1] "5 6 ExM.U vs ExN"
[1] 495   9
[1] 65
[1] "5 7 ExM.U vs InCGE"
[1] 506   9
[1] 66
[1] "5 8 ExM.U vs InMGE"
[1] 519   9
[1] 69
[1] "5 9 ExM.U vs IP"
[1] 537   9
[1] 70
[1] "5 10 ExM.U vs Mic"
[1] 567   9
[1] 72
[1] "5 11 ExM.U vs OPC"
[1] 569   9
[1] 72
[1] "5 12 ExM.U vs oRG"
[1] 590   9
[1] 72
[1] "5 13 ExM.U vs Per"
[1] 608   9
[1] 74
[1] "5 14 ExM.U vs PgG2M"
[1] 626   9
[1] 74
[1] "5 15 ExM.U vs PgS"
[1] 641   9
[1] 74
[1] "5 16 ExM.U vs vRG"
[1] 661   9
[1] 75
[1] "6 7 ExN vs InCGE"
[1] 665   9
[1] 76
[1] "6 8 ExN vs InMGE"
[1] 665   9
[1] 76
[1] "6 9 ExN vs IP"
[1] 667   9
[1] 76
[1] "6 10 ExN vs Mic"
[1] 676   9
[1] 76
[1] "6 11 ExN vs OPC"
[1] 676   9
[1] 76
[1] "6 12 ExN vs oRG"
[1] 690   9
[1] 77
[1] "6 13 ExN vs Per"
[1] 699   9
[1] 78
[1] "6 14 ExN vs PgG2M"
[1] 718   9
[1] 83
[1] "6 15 ExN vs PgS"
[1] 743   9
[1] 88
[1] "6 16 ExN vs vRG"
[1] 767   9
[1] 92
[1] "7 8 InCGE vs InMGE"
[1] 767   9
[1] 92
[1] "7 9 InCGE vs IP"
[1] 769   9
[1] 92
[1] "7 10 InCGE vs Mic"
[1] 773   9
[1] 92
[1] "7 11 InCGE vs OPC"
[1] 775   9
[1] 92
[1] "7 12 InCGE vs oRG"
[1] 795   9
[1] 97
[1] "7 13 InCGE vs Per"
[1] 808   9
[1] 100
[1] "7 14 InCGE vs PgG2M"
[1] 842   9
[1] 103
[1] "7 15 InCGE vs PgS"
[1] 879   9
[1] 104
[1] "7 16 InCGE vs vRG"
[1] 901   9
[1] 106
[1] "8 9 InMGE vs IP"
[1] 903   9
[1] 106
[1] "8 10 InMGE vs Mic"
[1] 907   9
[1] 106
[1] "8 11 InMGE vs OPC"
[1] 911   9
[1] 107
[1] "8 12 InMGE vs oRG"
[1] 947   9
[1] 111
[1] "8 13 InMGE vs Per"
[1] 961   9
[1] 114
[1] "8 14 InMGE vs PgG2M"
[1] 987   9
[1] 115
[1] "8 15 InMGE vs PgS"
[1] 1022    9
[1] 120
[1] "8 16 InMGE vs vRG"
[1] 1043    9
[1] 121
[1] "9 10 IP vs Mic"
[1] 1060    9
[1] 122
[1] "9 11 IP vs OPC"
[1] 1064    9
[1] 122
[1] "9 12 IP vs oRG"
[1] 1070    9
[1] 123
[1] "9 13 IP vs Per"
[1] 1081    9
[1] 123
[1] "9 14 IP vs PgG2M"
[1] 1085    9
[1] 124
[1] "9 15 IP vs PgS"
[1] 1094    9
[1] 126
[1] "9 16 IP vs vRG"
[1] 1098    9
[1] 126
[1] "10 11 Mic vs OPC"
[1] 1101    9
[1] 126
[1] "10 12 Mic vs oRG"
[1] 1125    9
[1] 127
[1] "10 13 Mic vs Per"
[1] 1125    9
[1] 127
[1] "10 14 Mic vs PgG2M"
[1] 1157    9
[1] 130
[1] "10 15 Mic vs PgS"
[1] 1198    9
[1] 136
[1] "10 16 Mic vs vRG"
[1] 1239    9
[1] 139
[1] "11 12 OPC vs oRG"
[1] 1243    9
[1] 141
[1] "11 13 OPC vs Per"
[1] 1243    9
[1] 141
[1] "11 14 OPC vs PgG2M"
[1] 1245    9
[1] 141
[1] "11 15 OPC vs PgS"
[1] 1252    9
[1] 141
[1] "11 16 OPC vs vRG"
[1] 1258    9
[1] 141
[1] "12 13 oRG vs Per"
[1] 1266    9
[1] 142
[1] "12 14 oRG vs PgG2M"
[1] 1268    9
[1] 142
[1] "12 15 oRG vs PgS"
[1] 1270    9
[1] 142
[1] "12 16 oRG vs vRG"
[1] 1270    9
[1] 142
[1] "13 14 Per vs PgG2M"
[1] 1288    9
[1] 145
[1] "13 15 Per vs PgS"
[1] 1305    9
[1] 146
[1] "13 16 Per vs vRG"
[1] 1313    9
[1] 146
[1] "14 15 PgG2M vs PgS"
[1] 1315    9
[1] 146
[1] "14 16 PgG2M vs vRG"
[1] 1318    9
[1] 146
[1] "15 16 PgS vs vRG"
[1] 1320    9
[1] 147
# Top hits p_val_adj < 0.05

length(unique(totaladj$geneId))
[1] 7
print(unique(totaladj$geneId))
[1] "CADM1" "GPM6A" "MAP1B" "APP"   "RTN4"  "DCTN2" "PFN2" 
write.table(totaladj, file="output/tables/TableS5D.txt", sep="\t", quote = F)
# Lower hits p_val < 0.05

#length(unique(total$geneId))
#print(unique(total$geneId))
#write.table(total, file="isoswitch.v2.txt", sep="\t", quote = F)