diff --git a/common/common/lin.py~ b/common/common/lin.py~ deleted file mode 100644 index 3cd910b4b85dba6e884fc51bbdd0c67b86643f41..0000000000000000000000000000000000000000 --- a/common/common/lin.py~ +++ /dev/null @@ -1,546 +0,0 @@ -from __future__ import print_function -from __future__ import division -from collections import defaultdict -from os.path import join -import os -import _pickle as cPickle -from seg import * -import analyseExprs as a -from itertools import repeat -import boolFs as b -import L1L2_cells_ids as layers - -linDataLoc = "../data/tracking/" - -def getSurvivingBetween(tss, linss, t1, t2): - linsL = filterLinssBetween(linss, (t1, t2)) - - cids = list() - - - for c in tss[t1]: - if tSearchN(linsL, c.cid, len(linsL), 0): - cids.append(c.cid) - - return cids - -def filterSurvivingUntil(tss, linss, t2): - for t, ts in tss.items(): - cids = getSurvivingBetween(tss, linss, t, t2) - filterT(ts, cids) - - return - -def fst(x): return x[0] - -def snd(x): return x[1] - -def filterT(ts, cids): - scids = set(cids) - ts.cells = dict([(c.cid, c) for c in ts if c.cid in scids]) - - neighsL1 = dict() - for cid, ds in ts.neighs.items(): - if cid in scids: - neighsL1[cid] = dict([(n, d) for n, d in ds.items() if n in scids]) - - ts.neighs = neighsL1 - -def filterL1(tss): - L1_cids = {10: layers.L1_10h, - 40: layers.L1_40h, - 96: layers.L1_96h, - 120: layers.L1_120h, - 132: layers.L1_132h} - - for t, ts in tss.items(): - filterT(ts, L1_cids[t]) - - return - - -def filterL1L2(tss): - L1_cids = {10: layers.L1_10h, - 40: layers.L1_40h, - 96: layers.L1_96h, - 120: layers.L1_120h, - 132: layers.L1_132h} - - L2_cids = {10: layers.L2_10h, - 40: layers.L2_40h, - 96: layers.L2_96h, - 120: layers.L2_120h, - 132: layers.L2_132h} - - for t, ts in tss.items(): - #scids1 = set(L1_cids[t]) - #scids2 = set(L2_cids[t]) - - filterT(ts, L1_cids[t] + L2_cids[t]) - -# ts.cells = dict([(c.cid, c) for c in ts -# if c.cid in scids1 or c.cid in scids2]) - - return - - -def getTs(linFn): - fn, _ = os.path.splitext(linFn) - (_, t1, _, t2) = fn.split("_") - - return (int(t1[:-1]), int(t2[:-1])) - - -def readLin(linFn): - with open(linFn, 'rb') as linF: - lin = cPickle.load(linF, encoding='latin1') - - return lin, getTs(linFn) - - -def readLins(linDataLoc): - lins = [readLin(join(linDataLoc, f)) for f in os.listdir(linDataLoc) - if os.path.isfile(join(linDataLoc, f))] - - return lins - - -def getLinMap(lin): - cMap = defaultdict(list) - for c1, c2 in lin[0]: - cMap[c2].append(c1) - - return cMap - - -def mkLinsSeries(lins): - return [(ts, getLinMap(lin)) - for lin, ts in sorted(lins, key=lambda lt: lt[1][0])] - - -def after(t, lt): - ts, _ = lt - (t1, t2) = ts - - return t >= t1 - - -def before(t, lt): - ts, _ = lt - (t1, t2) = ts - - return t <= t2 - - -def between(linSeries, t1, t2): - afters = [i for i, lt in enumerate(linSeries) if after(t1, lt)] - befores = [i for i, lt in enumerate(linSeries) if before(t2, lt)] - - s = afters[-1] - e = befores[0] - - return [(ts, lin) for ts, lin in linSeries[s:(e+1)]] - - -def mergeLins(lin1, lin2): - mLin = dict() - - for m, ds in lin1.items(): - mLin[m] = sum([lin2[d] for d in ds], []) - - return mLin - - -def mergeLinss1(lins, acc): - if not lins: return acc - - return mergeLinss1(lins[1:], mergeLins(acc, lins[0])) - - -def mergeLinss(lins): - if len(lins) > 2: return mergeLinss1(lins[2:], mergeLins(lins[0], lins[1])) - elif len(lins) == 2: return mergeLins(lins[0], lins[1]) - elif len(lins) == 1: return lins[0] - else: return [] - - -###################################################### -def concat(xss): - return sum(xss, []) - - -def tSearchN(mss, r, n, c): - if c == n: return [r] - else: - return concat([tSearchN(mss, d, n, c+1) for d in mss[c].get(r, [])]) - - -def tSearch1(mss, r): - try: - return mss[0][r] - except KeyError: - return [] - - -def tSearch2(mss, r): - try: - return sum([mss[1][d] for d in mss[0][r]], []) - except KeyError: - print(r) - return [] - - -def cvol(c): - if c: return c.vol - else: return 0.0 - - -def sumVol(cs, f=cvol): - return sum([cvol(c) for c in cs]) - - -def mergeCellVolsRefN(tss, mss): - volss = [] - - for i, ts in enumerate(tss): - volss.append(dict()) - for c in tss[0]: - volss[i][c.cid] = sum(cvol(ts.cells.get(ci, None)) for ci in tSearchN(mss, c.cid, i, 0)) - - return volss - - -def const1(c): - if c: return 1 - else: return 0 - - -def mergeCellNsRefN(tss, linss): - attrss = [] - - for i, ts in enumerate(tss): - attrss.append(dict()) - for c in tss[0]: - attrss[i][c.cid] = sum(const1(ts.cells.get(ci, None)) for ci in tSearchN(linss, c.cid, i, 0)) - - return attrss - - -def mergeCellVolsRef3(tss, mss): - (ts1, ts2, ts3) = tss - - vols1 = dict() - vols2 = dict() - vols3 = dict() - - for c in ts1: - vols1[c.cid] = c.vol - vols2[c.cid] = sumVol([ts2.cells[i] for i in tSearch1(mss, c.cid)]) - vols3[c.cid] = sumVol([ts3.cells[i] for i in tSearch2(mss, c.cid)]) - - return vols1, vols2, vols3 - - -def mergeCellVolsRef2(tss, mss): - (ts1, ts2) = tss - - vols1 = dict() - vols2 = dict() - - for c in ts1: - print(c.cid) - vols1[c.cid] = c.vol - vols2[c.cid] = sumVol([ts2.cells[i] for i in tSearch1(mss, c.cid)]) - - return vols1, vols2 - - -def divOr0(n1, n2): - if n2 == 0: return 0 - else: return n1 / n2 - - -def growthRates3(volss, m12): - vols1, vols2, vols3 = volss - gr2 = dict() - for i, vol in vols1.items(): - gr12 = divOr0((vols2[i] - vols1[i]), vols1[i]) - gr23 = divOr0((vols3[i] - vols2[i]), vols2[i]) - gr2[i] = 0.5 * (gr12 + gr23) - - gr22 = [zip(m12.get(i, []), repeat(gr, len(m12.get(i, [])))) - for i, gr in gr2.items()] - - return dict(sum(gr22, [])) - - -def growthRatesL2(volss): - vols1, vols2 = volss - gr1 = dict() - for i, vol in vols1.items(): - gr1[i] = (vols2[i] - vols1[i]) / vols1[i] - - return gr1 - - -def growthRatesR2(volss, m12): - vols1, vols2 = volss - gr2 = dict() - for i, vol in vols1.items(): - gr2[i] = (vols2[i] - vols1[i]) / vols1[i] - - gr22 = [zip(m12.get(i, []), repeat(gr, len(m12.get(i, [])))) - for i, gr in gr2.items()] - - return dict(sum(gr22, [])) - - -def mkSeries(): - d = "../data/organism" - dExprs = "../data/geneExpression/new_patterns" - - timepoints = [10, 40, 96, 120, 132] - lins = mkLinsSeries(readLins(linDataLoc)) - - tss = dict() - linss = dict() - - for t1, t2 in zip(timepoints, timepoints[1:]): - tss[t1] = a.readDataT1(d, dExprs, t1) - tss[t2] = a.readDataT1(d, dExprs, t2) - - linss[(t1, t2)] = mergeLinss(list(map(lambda x: x[1], between(lins, t1, t2)))) - - return tss, linss - -def mkSeriesIm0(): - d = "../data/yr01" - linDataLoc = "../data/yr01/tracking" - - timepoints = [0, 10, 40, 96, 120, 132] - lins = mkLinsSeries(readLins(linDataLoc)) - - tss = dict() - linss = dict() - - for t in timepoints: - tss[t] = a.readDataT1Im(d, t) - - for t1, t2 in zip(timepoints, timepoints[1:]): - linss[(t1, t2)] = mergeLinss(map(lambda x: x[1], - between(lins, t1, t2))) - - return tss, linss - - -def mkGrowthRates(): - tss, linss = mkSeries() - grs = dict() - - grs[10] = growthRatesL2(mergeCellVolsRef2((tss[10], tss[40]), - [linss[(10, 40)]])) - grs[40] = growthRates3(mergeCellVolsRef3((tss[10], tss[40], tss[96]), - [linss[(10, 40)], linss[(40, 96)]]), - linss[(10, 40)]) - grs[96] = growthRates3(mergeCellVolsRef3((tss[40], tss[96], tss[120]), - [linss[(40, 96)], linss[(96, 120)]]), - linss[(40, 96)]) - grs[120] = growthRates3(mergeCellVolsRef3((tss[96], tss[120], tss[132]), - [linss[(96, 120)], linss[(120, 132)]]), - linss[(96, 120)]) - grs[132] = growthRatesR2(mergeCellVolsRef2((tss[120], tss[132]), - [linss[(120, 132)]]), - linss[(120, 132)]) - - return grs - -def readGeomImT(fid, t): - ressFn = "../data/yr" + fid + "/resolutions.txt" - segFn = ("../data/yr" + - fid + - "/segmentations/YR_" + - fid + "_" + str(t) + - "h_segmented.tif") - - ress = readRess(ressFn) - img = readImage(segFn) - res = ress[t] - - ts = Tissue.fromImage(data=img, res=res) - for c in ts: - c.swapZX() - - return ts - -def mkSeriesGeomIm(fid): - ressFn = "../data/yr" + fid + "/resolutions.txt" - linDataLoc = "../data/yr" + fid + "/trackingFiles/" - ress = readRess(ressFn) - tpoints = sorted(ress.keys()) - - lins = mkLinsSeries(readLins(linDataLoc)) - tss = dict() - linss = dict() - - for t in tpoints: - tss[t] = readGeomImT(fid, t) - - for t1, t2 in zip(tpoints, tpoints[1:]): - linss[(t1, t2)] = mergeLinss(map(lambda x: - x[1], between(lins, t1, t2))) - - return tss, linss - - -################################### -def invertMap(ln): - iLin = [list(zip(ds, repeat(m, len(ds)))) for m, ds in ln.items()] - - return dict(sum(iLin, [])) - -def updateCells(ts, ts1, iLn): - for c1 in ts1: - try: - c = ts.cells[iLn[c1.cid]] - except KeyError: - print(c1.cid) - c = Cell(0, Vec3(0, 0, 0), 0, dict([(g, False) for g in geneNms])) - c1.exprs = dict([(g, c.exprs[g]) for g, v in c1.exprs.items()]) - - return - -def prev(t): - if t == 40: return 10 - elif t == 96: return 40 - elif t == 120: return 96 - elif t == 132: return 120 - else: return None - -def succ(t): - if t==10: return 40 - elif t==40: return 96 - elif t==96: return 120 - elif t==120: return 132 - else: return None - -def updateTs(tss, linss, tpoints): - #evolve patterns based on lineage - tssLin = {} - for t in tpoints[1:]: - tssLin[(prev(t), t)] = tss[t].updateExprsTs(tss[prev(t)], - invertMap(linss[(prev(t), t)])) - - return tssLin - -def updateTsBF(tss, fs): - return dict((t, tss[t].updateExprsBF(fs[t])) for t in tss.keys()) - -def toNewmanFn(fn): - (fn, ext) = os.path.splitext(fn) - return fn + "000" + ".tif" - -def plotNewman1(tss, i, geneNm, outDir="."): - from subprocess import call - import os - - confFile = "newman.conf" - visCmd = "/home2/Organism/tools/plot/bin/newman" - mergeCmd = "convert" - - b.cellsToDat("ts1.dat", list(tss)) - - call([visCmd, - "-shape", "sphereVolume", - "-d", "3", - "ts1.dat", - "-column", str(i+1), - "-output", "tiff", - geneNm, - "-camera", confFile, - "-min", str(0), - "-max", str(1)]) - - os.remove("ts1.dat") - - return - -def plotNewman2(tss, tss1, i, geneNm, outDir="."): - from os.path import join - from subprocess import call - import os - - confFile = "newman.conf" - visCmd = "/home2/Organism/tools/plot/bin/newman" - mergeCmd = "convert" - - b.cellsToDat("ts1.dat", list(tss)) - b.cellsToDat("ts2.dat", list(tss1)) - - call([visCmd, - "-shape", "sphereVolume", - "-d", "3", - "ts1.dat", - "-column", str(i+1), - "-output", "tiff", - "ts1", - "-camera", confFile, - "-min", str(0), - "-max", str(1)]) - - call([visCmd, - "-shape", "sphereVolume", - "-d", "3", - "ts2.dat", - "-column", str(i+1), - "-output", "tiff", - "ts2", - "-camera", confFile, - "-min", str(0), - "-max", str(1)]) - - call([mergeCmd, - toNewmanFn("ts1.tif"), - toNewmanFn("ts2.tif"), - "+append", - "-pointsize", "20", - "-fill", "white", - "-annotate", "+10+40", geneNm, - join(outDir, geneNm+".tif")]) - - os.remove("ts1.dat") - os.remove("ts2.dat") - os.remove(toNewmanFn("ts1.tif")) - os.remove(toNewmanFn("ts2.tif")) - - -def tFn(t, ref=""): - return "ts_t" + str(t) + ref + ".txt" - - -def toOrgs(tss, ref="", outDir="."): - for t, ts in tss.items(): - b.cellsToOrg(join(outDir, tFn(t, ref)), list(ts)) - - return - -def snd(x): - return x[1] - -def betweenTs(t, bounds): - ts, te = bounds - return t >= ts and t <= te - -def filterLinssBetween(linss, bounds): - t1, t2 = bounds - linssL = sorted([(ts, ln) for ts, ln in linss.items() - if betweenTs(ts[0], (t1, t2)) and betweenTs(ts[1], (t1, t2))], key=lambda p: p[0][0]) - - return list(map(snd, linssL)) - -def filterTssBetween(tss, tbounds): - t1, t2 = tbounds - - tssL = sorted([(t, ts) for t, ts in tss.items() if betweenTs(t, (t1, t2))], key=fst) - - return list(map(snd, tssL)) - diff --git a/common/common/seg.py~ b/common/common/seg.py~ deleted file mode 100644 index 8fe48542254b939da38243fc05e7587c0dc0170e..0000000000000000000000000000000000000000 --- a/common/common/seg.py~ +++ /dev/null @@ -1,601 +0,0 @@ -from __future__ import division -from __future__ import print_function -import numpy as np -import ast -from tifffile import TiffFile -import matplotlib.pyplot as plt -import boolFs as b -from copy import deepcopy -from mpl_toolkits.mplot3d import axes3d - - -geneNms = ['AG', - 'AHP6', - 'ANT', - 'AP1', - 'AP2', - 'AP3', - 'AS1', - 'ATML1', - 'CLV3', - 'CUC1_2_3', - 'ETTIN', - 'FIL', - 'LFY', - 'MP', - 'PHB_PHV', - 'PI', - 'PUCHI', - 'REV', - 'SEP1', - 'SEP2', - 'SEP3', - 'STM', - 'SUP', - 'SVP', - 'WUS'] - -class Vec3(object): - def __init__(self, x, y, z): - self.x = x - self.y = y - self.z = z - - def dist(self, v1): - return np.sqrt((self.x - v1.x)**2 + - (self.y - v1.y)**2 + - (self.z - v1.z)**2) - - def toOrganism(self, sep=' '): - return sep.join([str(self.x), str(self.y), str(self.z)]) - - def swapZX(self): - xt = self.x - self.x = self.z - self.z = xt - - def toArray(self): - return np.array([self.x, self.y, self.z]) - - def __repr__(self): - return self.toOrganism() - -class Cell(object): - def __init__(self, cid, pos, vol, exprs=None): - self.cid = cid #Int - self.pos = pos #Vec3 - self.vol = vol #Real - self.exprs = exprs #Map String Bool - - def dist(self, c1): - return self.pos.dist(c1.pos) - - def addExprs(self, exprs, geneNms): - for gn in geneNms: - if not gn in exprs.keys(): - exprs[gn] = False - - self.exprs = exprs - self.geneNms = geneNms - - return - - def updateExprs(self, gfs): - gexprsUpd = dict() - - for g, _ in self.exprs.items(): - if g in gfs: - gf = gfs[g] - (acts, reprs, out, fsAct, fsRepr) = gf - inputs = acts + reprs - vals = [self.exprs[i] for i in inputs] - gexprsUpd[g] = b.evalI((acts, reprs, out), vals, fsAct, fsRepr) - else: - gexprsUpd[g] = self.exprs[g] - - ncell = Cell(self.cid, self.pos, self.vol) - ncell.addExprs(gexprsUpd, deepcopy(self.geneNms)) - - return ncell - - def updateExprsC(self, c): - exprsUpd = dict([(g, c.exprs[g]) for g, v in self.exprs.items()]) - nc = Cell(self.cid, self.pos, self.vol, exprsUpd) - nc.geneNms = self.geneNms - return nc - - def _boolToInt(self, b): - if b: return 1 - - return 0 - - def _exprsToOrg(self, sep=' '): - exprs = [str(self._boolToInt(self.exprs[gn])) for gn in self.geneNms] - - return sep.join(exprs) - - def swapZX(self): - self.pos.swapZX() - - @classmethod - def fromTV(cls, nelems): - def parsePos(spos): - pos = list(np.fromstring(spos.replace('[', '').replace(']', ''), - sep=' ')) - return pos - - cid = int(nelems['cid']) - vol = float(nelems['volume']) - [xc, yc, zc] = parsePos(nelems['center']) - pos = Vec3(xc, yc, zc) - - return cls(cid, pos, vol) - - @classmethod - def fromImage(cls, cd, res): - volVoxel = np.prod(res) - xc, yc, zc = cd.centroid - nVoxels = cd.area - - xc, yc, zc = np.multiply(np.array([xc, yc, zc]), res) - - vol = nVoxels * volVoxel - pos = Vec3(xc, yc, zc) - cid = cd.label - - return cls(cid, pos, vol) - - def toOrganism(self): - if self.exprs: - return (' '.join([self.pos.toOrganism(), str(self.vol), - str(self.cid), self._exprsToOrg()])) - else: - return (' '.join([self.pos.toOrganism(), str(self.vol), - str(self.cid)])) - - def toText(self, sep=" "): - if self.exprs: - return (sep.join([str(self.cid), - self.pos.toOrganism(sep=sep), - str(self.vol), - self._exprsToOrg(sep=sep)])) - else: - return sep.join([str(self.cid), self.pos.toOrganism(sep=sep)]) - - def centreToArray(self): - return self.pos.toArray() - - def getOnGenes(self): - return frozenset([g for g, st in self.exprs.items() if st]) - frozenset(['WUS']) - - def __repr__(self): - return self.toOrganism() - -class Tissue(object): - def __init__(self, cells, neighs, ecellIds = None): - self.cells = cells #Map Int Cell - self.neighs = neighs #Map Int (Map Int Dist) - self.L1_cells = ecellIds #Set Int - - def __iter__(self): - #iterate over a list instead of map - return iter([v for _, v in self.cells.items()]) - - @classmethod - def _extractCells(cls, data, res): - from skimage.measure import regionprops - cellsSeg = regionprops(data) - - cells = dict() - for cd in cellsSeg: - cells[int(cd.label)] = Cell.fromImage(cd, res) - - return cells - - @classmethod - def _getConn(cls, cells, data): - from skimage.future.graph import RAG - - connG = RAG(data) - neigh = dict() - - for n in connG.nodes: - nbrs = connG.adj[n].keys() - ds = [cells[n].dist(cells[nb]) for nb in nbrs] - - neigh[n] = dict(zip(nbrs, ds)) - - return neigh - - @classmethod - def fromTV(cls, fnGeom, fnNeigh): - def parseLn(ln): - ln = str.strip(ln) - nelems = dict() - elems = str.split(ln, ',') - for elm in elems: - [name, val] = str.split(elm, ':') - name = str.strip(name) - val = str.strip(val) - nelems[name] = val - - return nelems - - cells = dict() - neighs = dict() - ecellIds = [] - - with open(fnGeom, 'r') as f: - for ln in f: - nelems = parseLn(ln) - c = Cell.fromTV(nelems) - cells[c.cid] = c - - with open(fnNeigh, 'r') as f: - for ln in f: - nelems = ast.literal_eval(ln) - cid = nelems['cid'] - ns = nelems['neighbors ids'] - - if 1 in ns: - ecellIds.append(cid) - ns.remove(1) - - neighs[cid] = dict((n, cells[cid].dist(cells[n])) for n in ns) - - return cls(cells, neighs, set(ecellIds)) - - @classmethod - def fromImage(cls, data, res): - cells = cls._extractCells(data, res) - conn = cls._getConn(cells, data) - - return cls(cells, conn) - - def adjView(self): - return dict([(cid, ds.keys()) for cid, ds in self.neighs.items()]) - - def relabelFromMap(self, relabelM): - neighsR = [] - - for c in self: - c.cid = relabelM[c.cid] - - for c, dists in self.neighs.items(): - reNs = dict([(relabelM[n], d) for n, d in dists.items()]) - neighsR.append((relabelM[c], reNs)) - - self.neighs = dict(neighsR) - - return - - def relabel(self): - relabelM = dict([(c.cid, i) for i, c in enumerate(self)]) - - self.relabelFromMap(relabelM) - - return - - def filterL1(self): - if not self.neighs: - return - - L1_cell_ids = self.neighs[1] - - L1_cells = dict([(c.cid, c) for c in self if c.cid in L1_cell_ids]) - L1_neighs = dict([(cid, nbrs) for cid, nbrs in deepcopy(self.neighs.items()) - if cid in L1_cells_ids]) - - return Tissue(L1_cells, L1_neighs) - - def toNeigh(self): - - def cellToNeigh(cid, nbs): - nbsIds = nbs.keys() - nNbs = len(nbsIds) - nbsDists = [nbs[nid] for nid in nbsIds] - - return " ".join([str(cid), str(nNbs)] - + map(str, nbsIds) - + map(str, nbsDists)) - - ncells = sum(1 for _ in self) - header = str(ncells) + " " + "1" - - cellReprs = ([cellToNeigh(cid, nbs) - for cid, nbs in self.neighs.items()]) - - return "\n".join([header] + cellReprs) - - def toOrganism(self): - headerComment = "#x y z vol id" - ns = str(len(self.cells)) + " 5" - content = '\n'.join([cell.toOrganism() for cell in self]) - - return "\n".join([headerComment, ns, content]) - - def centresToArray(self): - return np.array([c.centreToArray() for c in self]) - - -class STissue(Tissue): - def __init__(self, cells, neighs, gexprs=None, geneNms=None): - super(STissue, self).__init__(cells, neighs) - if gexprs and geneNms: - self.addExprs(gexprs, geneNms) - - def _defGeneExprs(self): - return dict([(gn, False) for gn in self.geneNms]) - - def addExprs(self, gexprs, geneNms): - self.geneNms = geneNms - for cell in self: - try: - cell.addExprs(gexprs[cell.cid], geneNms) - except KeyError: - cell.addExprs(self._defGeneExprs(), geneNms) - - return - - def relabel(self): - super(STissue, self).relabel() - - @classmethod - def fromImage(cls, data, res, gexprs, geneNms): - ts = Tissue.fromImage(data, res) - return cls(ts.cells, ts.neighs, gexprs, geneNms) - - @classmethod - def fromTV(cls, fnGeom, fnNeigh, gexprs, geneNms): - ts = Tissue.fromTV(fnGeom, fnNeigh) - return cls(ts.cells, ts.neighs, gexprs, geneNms) - - def toOrganism(self): - headerComment = "#x y z vol id " + ' '.join(self.geneNms) - ns = str(len(self.cells)) + " " + str(len(self.geneNms)+5) - content = '\n'.join([cell.toOrganism() for cell in self]) - - return "\n".join([headerComment, ns, content]) - - def updateExprsBF(self, gfs): - ncells = dict([(c.cid, c.updateExprs(gfs)) for c in self]) - return STissue(ncells, deepcopy(self.neighs)) - - def updateExprsTs(self, ts, iLn): - ncells = {} - for c1 in self: - try: - c = ts.cells[iLn[c1.cid]] - except KeyError: - c = Cell(0, Vec3(0, 0, 0), 0, dict([(g, False) for g in geneNms])) - - ncells[c1.cid] = c1.updateExprsC(c) - - return STissue(ncells, deepcopy(self.neighs)) - - def toDict(self): - tsDict = {} - - for c in self: - tsDict[c.cid] = dict(((g, b.boolToInt(v)) for g, v in c.exprs)) - - return tsDict - - def toCSV(self): - headerComment = "id,x,y,z,vol," + ','.join(self.geneNms) - - content = '\n'.join([cell.toText(sep=",") for cell in self]) - - return "\n".join([headerComment, content]) - - - - -#plotting -class IndexTracker(object): - def __init__(self, ax, X): - self.ax = ax - ax.set_title('use scroll wheel to navigate images') - - self.X = X - self.slices, rows, cols = X.shape - self.ind = self.slices//2 - - self.im = ax.imshow(self.X[self.ind, :, :]) - self.update() - - def onscroll(self, event): - if event.button == 'up': - self.ind = (self.ind + 1) % self.slices - else: - self.ind = (self.ind - 1) % self.slices - self.update() - - def update(self): - self.im.set_data(self.X[self.ind, :, :]) - self.ax.set_ylabel('slice %s' % self.ind) - self.im.axes.figure.canvas.draw() - -def plot3d(data): - fig, ax = plt.subplots(1, 1) - tracker = IndexTracker(ax, data) - fig.canvas.mpl_connect('scroll_event', tracker.onscroll) - plt.show() - - -############### reading -def getResolution(f): - try: - x1, x2 = f.pages.pages[0].tags['XResolution'] - y1, y2 = f.pages.pages[0].tags['YResolution'] - xres = x2 / x1 - yres = y2 / y1 - zres = 1.0 - except: - xres = 1.0 - yres = 1.0 - zres = 1.0 - - return np.array([xres, yres, zres]) - -def readImage(fseg): - f = TiffFile(fseg) - seg = f.asarray().astype(int) - - return seg - -def readImageRegionPoints(fseg, res=[1, 1, 1], d=2): - from skimage.measure import regionprops - - cell_regions = dict() - - seg = readImage(fseg) - cellsSeg = regionprops(seg) - - for cregion in cellsSeg: - cell_regions[cregion.label] = np.multiply(cregion.coords[:, :d], res[:d]) - - return cell_regions - -def readRegionVolume(fseg, res=[1, 1, 1]): - from skimage.measure import regionprops - - cell_regions = dict() - - seg = readImage(fseg) - cellsSeg = regionprops(seg) - - for cregion in cellsSeg: - cell_regions[cregion.label] = cregion.area*np.product(res) - - return cell_regions - - - -def mergePoints(cell_regions, cids): - return np.vstack(tuple([cell_regions[cid] for cid in cids])) - -def getCellRegions(fid="01", d=2): - segFn = "../data/yr{fid}/segmentations/YR_{fid}_{t}h_segmented.tif" - ressFn = "../data/yr" + fid + "/resolutions.txt" - ress = readRess(ressFn) - tpoints = sorted(ress.keys()) - - fsegs = dict([(t, segFn.format(fid=fid, t=str(t))) for t in tpoints]) - - return dict([(t, readImageRegionPoints(fn, res=ress[t], d=d)) - for t, fn in fsegs.items()]) - -def getCellRegionVolumes(fid="01"): - fsegs = {10: "../data/yr01/segmentations/YR_01_10h_segmented.tif", - 40: "../data/yr01/segmentations/YR_01_40h_segmented.tif", - 96: "../data/yr01/segmentations/YR_01_96h_segmented.tif", - 120: "../data/yr01/segmentations/YR_01_120h_segmented.tif", - 132: "../data/yr01/segmentations/YR_01_132h_segmented.tif"} - - return dict([(t, readRegionVolume(fn, res=ress[t])) - for t, fn in fsegs.items()]) - -def readExprs(fn): - #returns: Map Int (Map String Bool) Cid -> (GeneName -> {On/Off}) - gexprs = dict() - - with open(fn) as f: - geneNms = f.readline().rstrip().split(' ')[1:] - - for line in f: - els = line.rstrip().split(' ') - cid = els[0] - exprs = els[1:] - gexprs[int(cid)] = dict(zip(geneNms, [bool(int(exp)) - for exp in exprs])) - - return gexprs - -def writeOrg(fnGeom, fnNeigh, ts): - with open(fnGeom, 'w+') as f: - f.write(ts.toOrganism()) - - with open(fnNeigh, 'w+') as f: - f.write(ts.toNeigh()) - - return - -def readRess(fn): - import ast - - with open(fn, 'r') as f: - sress = ("{" + f.read().replace("=", ":").replace("h", ""). - replace("\n", "").replace("t", "") + "}") - ress = ast.literal_eval(sress) - - return ress - -#functions for common tasks -def mkOrgs(d, dExprs, dOut): - from os.path import join - - timepoints = [10, 40, 96, 120, 132] - geomF = "_segmented_tvformat_0_volume_position.txt" - neighF = "_segmented_tvformat_0_neighbors_py.txt" - exprF = "h.txt" - - initF = ".init" - orgNeighF = ".neigh" - - for t in timepoints: - print(t) - fnGeom = join(d, "t" + str(t) + geomF) - fnNeigh = join(d, "t" + str(t) + neighF) - fnExpr = join(dExprs, "t_" + str(t) + exprF) - fnGeomOut = join(dOut, "t" + str(t) + initF) - fnNeighOut = join(dOut, "t" + str(t) + orgNeighF) - - gexprs = readExprs(fnExpr) - ts = STissue.fromTV(fnGeom, fnNeigh, gexprs, geneNms) - - writeOrg(fnGeomOut, fnNeighOut, ts) - - return - -def mkOrgGeomState(fim, fseg, fexpr, fout): - (data, res) = readImage(fim, fseg) - gexprs = readExprs(fexpr) - - ts = STissue.fromImage(data, res, gexprs, geneNms) - - writeOrg(fout, ts) - -def mkOrgGeom(fim, fseg, fout): - (data, res) = readImage(fim, fseg) - ts = Tissue.fromImage(data, res) - - writeOrg(fout, ts) - - -def plotTissue(ts, hl): - xs = [c.pos.x for c in ts] - ys = [c.pos.y for c in ts] - zs = [c.pos.z for c in ts] - - cs = ['red' if c.cid in hl else 'blue' for c in ts] - - fig = plt.figure() - ax = fig.add_subplot(211, projection='3d') - ax.view_init(elev=95, azim=270) - ax.scatter(xs, ys, zs, c=cs, alpha=0.9, s=200, - edgecolors='black', linewidths=0.5) - - plt.show() - -def filterL1(ts): - if not ts.neighs: - return - - L1_cell_ids = ts.neighs[1] - - L1_cells = dict([(c.cid, c) for c in ts if c.cid in L1_cell_ids]) - L1_neighs = dict([(cid, nbrs) for cid, nbrs in deepcopy(list(ts.neighs.items())) - if cid in L1_cell_ids]) - - return Tissue(L1_cells, L1_neighs) diff --git a/common/setup.py~ b/common/setup.py~ deleted file mode 100644 index e69de29bb2d1d6434b8b29ae775ad8c2e48c5391..0000000000000000000000000000000000000000