I have matrix file which is basically a spearman correlation matrix between genes across various cell type. So now Im trying to find out which set of genes or group of genes whose correlation value is lets say greater than 0.6 if I set that as my threshold. How can I do that? I'm posting a subset of my data. It's a 502 x 502 matrix.
ACTL6B ACTR5 ACTR6
ACTL6B 1 0.6 -0.4
ACTR5 0.4 1 -0.3
ACTR6 -0.4 -0.3 1
So I don't want correlation between same set of genes which would be 1. I want other comparison. Like, lets say, ACTL6B and ACTR5 whose correlation is 0.6. I would like to keep those values and genes.
Any help or suggestion would be highly appreciated
Nice idea, though you would need 2 as margin, else the results will be transposed (as you see in your example)
as.data.frame(apply(test, 2, function(x) ifelse (abs(x) >=0.6,x,"NA")))
A simpler solution is to just use sapply
as.data.frame(sapply(test, function(x) ifelse (abs(x) >=0.6,x,"NA")))
sapply doesn't save/print row names.
It does!
You could use something like Python for this:
#!/usr/bin/env python
import sys
ctr = 0
genes = []
map = {}
for line in sys.stdin:
elems = [x for x in line.strip().split()]
if ctr == 0:
genes = elems[1:]
for gene in genes:
map[gene] = {}
else:
paired_gene = elems[0]
scores = [float(x) for x in elems[1:]]
for idx, score in enumerate(scores):
gene = genes[idx]
map[gene][paired_gene] = score
map[paired_gene][gene] = score
ctr += 1
for gene in genes:
for paired_gene in genes:
if gene != paired_gene and map[gene][paired_gene] >= 0.6:
sys.stdout.write("%s\t%s\t%f\n" % (gene, paired_gene, map[gene][paired_gene]))
d = data.frame(matrix(1:15, nrow=3))
# linearize the data
l = unlist(d)
# Put whatever filters you like
l[l>10] = NA
# re-assmeble
d = data.frame(matrix(l, nrow = dim(d)[1]))
This is the neatest!!!
Yes then you can deal with them as you want, e.g. many functions have an na.rm or na.omit option. Actually, this solution is wrong because it has "NA". The quotes makes it a string of value NA. What you want is NA (no quote) which is of logical type, not a string. Note that I gave you this solution right from the start.
Log in to answer this question.
This is not a bioinformatics question but rather an R programming question. Anyway try something along the line of
but next time post on StackOverflow.
Solution provided above by Jean-Karim Heriche is better and shorter :)
+1. Sorry, I missed it in the first instance.