-
Notifications
You must be signed in to change notification settings - Fork 10
/
sc_genestat.m
55 lines (45 loc) · 1.21 KB
/
sc_genestat.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
function [lgu, dropr, lgcv, genes, X, ...
removedgeneidx, removedT] = sc_genestat(X, genelist, ...
sortit, removeinf)
if nargin < 4, removeinf = true; end
if nargin < 3, sortit = true; end
if nargin < 2 || isempty(genelist)
genelist = "Gene"+string(1:size(X,1)).';
end
genelist=genelist(:);
geneidx = 1:length(genelist);
dropr = 1 - sum(X > 0, 2) ./ size(X, 2);
% m = X./sum(X);
% lgm = log(std(m, [], 2, 'omitnan') ./ mean(m, 2, 'omitnan'));
u = mean(X, 2, 'omitnan');
cv = std(X, [], 2, 'omitnan') ./ u;
lgu = log(u+1);
lgcv = log(cv+1);
genes = genelist;
if sortit
[xyz, si] = sortrows([lgu, dropr, lgcv], [1, 3, 2]);
lgu = xyz(:, 1);
dropr = xyz(:, 2);
lgcv = xyz(:, 3);
genes = genelist(si);
X=X(si,:);
geneidx = geneidx(si);
end
si = isnan(lgu) | isinf(lgu) | isnan(lgcv) | isinf(lgcv);
if removeinf && any(si)
removedT = table(genes(si), lgu(si), lgcv(si), dropr(si), ...
zeros(sum(si),1), ...
ones(sum(si),1), ...
ones(sum(si),1),...
zeros(sum(si),1));
lgu(si) = [];
lgcv(si) = [];
dropr(si) = [];
genes(si) = [];
X(si, :) = [];
removedgeneidx = geneidx(si);
else
removedgeneidx = [];
removedT = [];
end
end