Computes RootedGraph.expr(), RootedGraph.delta(), RootedGraph.connectivity(), RootedGraph.connectivity_pvalue(), RootedGraph.centroid() and Benjamini-Holchberg adjusted q-values for all genes that are expressed in more than 'filtercells' cells and whose maximum expression
(self, n=500, filtercells=0, filterexp=0.0, annotation={})
| 1866 | return self.po |
| 1867 | |
| 1868 | def save(self, n=500, filtercells=0, filterexp=0.0, annotation={}): |
| 1869 | """ |
| 1870 | Computes RootedGraph.expr(), RootedGraph.delta(), RootedGraph.connectivity(), |
| 1871 | RootedGraph.connectivity_pvalue(), RootedGraph.centroid() and Benjamini-Holchberg adjusted q-values for all |
| 1872 | genes that are expressed in more than 'filtercells' cells and whose maximum expression |
| 1873 | value is above 'filterexp'. The optional argument 'annotation' allos to include a dictionary |
| 1874 | with lists of genes to be annotated in the table. The output is stored in a tab separated |
| 1875 | file called name.genes.txt. |
| 1876 | """ |
| 1877 | pol = [] |
| 1878 | with open(self.name + '.genes.tsv', 'w') as ggg: |
| 1879 | cul = 'Gene\tCells\tMean\tMin\tMax\tConnectivity\tp_value\tq-value (BH)\tCentroid\tDispersion\t' |
| 1880 | for m in sorted(annotation.keys()): |
| 1881 | cul += m + '\t' |
| 1882 | ggg.write(cul[:-1] + '\n') |
| 1883 | lp = sorted(self.dicgenes.keys()) |
| 1884 | for gi in lp: |
| 1885 | if self.expr(gi) > filtercells and self.delta(gi)[2] > filterexp: |
| 1886 | pol.append(self.connectivity_pvalue(gi, n=n)) |
| 1887 | por = benjamini_hochberg(pol) |
| 1888 | mj = 0 |
| 1889 | for gi in lp: |
| 1890 | po = self.expr(gi) |
| 1891 | m1, m2, m3 = self.delta(gi) |
| 1892 | p1, p2 = self.centroid(gi) |
| 1893 | if po > filtercells and m3 > filterexp: |
| 1894 | cul = gi + '\t' + str(po) + '\t' + str(m1) + '\t' + str(m2) + '\t' + str(m3) + '\t' +\ |
| 1895 | str(self.connectivity(gi)) + '\t' + str(pol[mj]) + '\t' + str(por[mj]) +\ |
| 1896 | '\t' + str(p1) + '\t' + str(p2) + '\t' |
| 1897 | for m in sorted(annotation.keys()): |
| 1898 | if gi in annotation[m]: |
| 1899 | cul += 'Y' + '\t' |
| 1900 | else: |
| 1901 | cul += 'N' + '\t' |
| 1902 | ggg.write(cul[:-1] + '\n') |
| 1903 | mj += 1 |
| 1904 | centr = [] |
| 1905 | disp = [] |
| 1906 | centr2 = [] |
| 1907 | disp2 = [] |
| 1908 | f = open(self.name + '.genes.tsv', 'r') |
| 1909 | for n, line in enumerate(f): |
| 1910 | if n > 0: |
| 1911 | sp = line[:-1].split('\t') |
| 1912 | if float(sp[7]) <= 0.05: |
| 1913 | centr.append(float(sp[1])) |
| 1914 | disp.append(float(sp[5])) |
| 1915 | else: |
| 1916 | centr2.append(float(sp[1])) |
| 1917 | disp2.append(float(sp[5])) |
| 1918 | f.close() |
| 1919 | pylab.scatter(centr2, disp2, alpha=0.2, s=9, c='b') |
| 1920 | pylab.scatter(centr, disp, alpha=0.3, s=9, c='r') |
| 1921 | pylab.xlabel('cells') |
| 1922 | pylab.ylabel('connectivity') |
| 1923 | pylab.yscale('log') |
| 1924 | pylab.ylim(0.01, 1) |
| 1925 | pylab.xlim(0, max(centr+centr2)) |
nothing calls this directly
no test coverage detected