1 #include "SoloFeature.h"
2 #include "streamFuns.h"
3 #include "Stats.h"
4 #include "GlobalVariables.h"
5 
statsOutput()6 void SoloFeature::statsOutput()
7 {
8     ofstream &strOut=ofstrOpen(outputPrefix+"Summary.csv", ERROR_OUT, P);
9     //Sequencing
10     strOut << "Number of Reads," << g_statsAll.readN <<'\n';
11     strOut << "Reads With Valid Barcodes," << 1.0 - double( readBarSum->stats.numInvalidBarcodes() + readFeatSum->stats.numInvalidBarcodes() )/g_statsAll.readN <<'\n';
12     strOut << "Sequencing Saturation," << readFeatSum->stats.numSequencingSaturation() <<'\n';
13 
14     if (pSolo.type != pSolo.SoloTypes::SmartSeq) {//qualHist for CB+UMI
15         uint64 q30=0, ntot=0;
16         for (uint32 ix=0; ix<256; ix++) {
17             ntot += readBarSum->qualHist[ix];
18             if (ix >= (P.readQualityScoreBase + 30))
19                 q30 += readBarSum->qualHist[ix];
20         };
21 
22         strOut << "Q30 Bases in CB+UMI,"   << double(q30)/ntot <<'\n';
23     };
24 
25     {//qualHist for RNA
26         uint64 q30=0, ntot=0;
27         for (int ichunk=0; ichunk<P.runThreadN; ichunk++) {
28             for (uint32 imate=0; imate<P.readNmates; imate++) {
29                 for (uint32 ix=0; ix<256; ix++) {
30                     ntot += RAchunk[ichunk]->RA->qualHist[imate][ix];
31                     if (ix >= (P.readQualityScoreBase + 30))
32                         q30 += RAchunk[ichunk]->RA->qualHist[imate][ix];
33                 };
34             };
35         };
36         strOut << "Q30 Bases in RNA read," << double(q30)/ntot <<'\n';
37     };
38 
39     strOut << "Reads Mapped to Genome: Unique+Multiple," << double(g_statsAll.mappedReadsU+g_statsAll.mappedReadsM)/g_statsAll.readN <<'\n';
40     strOut << "Reads Mapped to Genome: Unique," << double(g_statsAll.mappedReadsU)/g_statsAll.readN <<'\n';
41 
42     string mapfeat=SoloFeatureTypes::Names[featureType];
43     //if (featureType==SoloFeatureTypes::Gene || featureType==SoloFeatureTypes::GeneFull)
44     //    mapfeat="Transcriptome";
45 
46     strOut << "Reads Mapped to "<< mapfeat << ": Unique+Multipe " << SoloFeatureTypes::Names[featureType] <<"," << double( readFeatSum->stats.numMappedToTranscriptome() )/g_statsAll.readN <<'\n';
47     strOut << "Reads Mapped to "<< mapfeat << ": Unique " << SoloFeatureTypes::Names[featureType] <<"," << double( readFeatSum->stats.numMappedToTranscriptomeUnique() )/g_statsAll.readN <<'\n';
48 
49     if (pSolo.cellFilter.type[0]!="None" && (featureType==SoloFeatureTypes::Gene || featureType==SoloFeatureTypes::GeneFull)) {
50         //if (pSolo.cellFilter.type[0]=="CellRanger2.2")
51         {
52             strOut << "Estimated Number of Cells," << filteredCells.nCells <<'\n';
53 
54             strOut << "Unique Reads in Cells Mapped to " << SoloFeatureTypes::Names[featureType] << "," << filteredCells.nReadInCellsUnique <<'\n';
55             strOut << "Fraction of Unique Reads in Cells," << double(filteredCells.nReadInCellsUnique) / readFeatSum->stats.numMappedToTranscriptomeUnique() <<'\n';
56             strOut << "Mean Reads per Cell," << filteredCells.meanReadPerCellUnique <<'\n';
57             strOut << "Median Reads per Cell," << filteredCells.medianReadPerCellUnique <<'\n';
58 
59             strOut << "UMIs in Cells," << filteredCells.nUMIinCells <<'\n';
60             strOut << "Mean UMI per Cell," << filteredCells.meanUMIperCell <<'\n';
61             strOut << "Median UMI per Cell," << filteredCells.medianUMIperCell <<'\n';
62 
63             strOut << "Mean "   << SoloFeatureTypes::Names[featureType] << " per Cell," << filteredCells.meanGenePerCell <<'\n';
64             strOut << "Median " << SoloFeatureTypes::Names[featureType] << " per Cell," << filteredCells.medianGenePerCell <<'\n';
65             strOut << "Total "  << SoloFeatureTypes::Names[featureType] << " Detected," << filteredCells.nGeneDetected <<'\n';
66         };
67 
68         //output UMI per cell, sorted
69         ofstream &strOutUMIperCell = ofstrOpen(outputPrefix+"UMIperCellSorted.txt", ERROR_OUT, P);
70 
71         for (auto & n : nUMIperCBsorted) {
72             if (n==0)
73                 break;
74             strOutUMIperCell << n <<'\n';
75         };
76         strOutUMIperCell.close();
77     };
78 
79     strOut.close();
80 };