I have what I imagine is a very simple error regarding the use of WGCNA in R. I am trying to create a weighted correlation network by setting a hard threshold using the function pickHardThreshold. However, when I pass in my expression set I get a strange error:
Error in lm.fit(x, y, offset = offset, singular.ok = singular.ok, ...) :
NA/NaN/Inf in 'x'
I imagine there is something wrong with my understanding of how to use the function but there seems to be no tutorials or examples where a hard threshold is set (there are many based on soft thresholds, such as here).
I am passing in an expression set matrix with the samples as the rows and the genes as the columns as described in the documentation for the function
expression data in a matrix or data frame. Rows correspond to samples and columns to genes.
Also, there are no NA's in my dataset which was my first check based on the error. Just to give a reproducible example of the error see below:
exp_matrix<-matrix(data=runif(100),ncol=10,
dimnames=list(paste0("sample_",1:10),
paste0("gene_",1:10)))
pickHardThreshold(exp_matrix)
Error in lm.fit(x, y, offset = offset, singular.ok = singular.ok, ...) :
NA/NaN/Inf in 'x'
Any help would be greatly appreciated!
Thanks.
•
written
by
Al Murphy •
How big is the original dataset? I am asking because by increasing the number of genes, the function works just fine
> exp_matrix<-matrix(data=runif(10000),ncol=1000,
+ dimnames=list(paste0("sample_",1:10),
+ paste0("gene_",1:1000)))
> pickHardThreshold(exp_matrix)
Cut p.value SFT.R.sq slope. truncated.R.2 mean.k. median.k. max.k.
1 0.10 0.77000 0.129000 18.000 0.9750 781.000 781.0 825
2 0.15 0.66000 0.113000 -13.000 0.9570 676.000 676.0 736
3 0.20 0.55500 0.160000 -11.600 0.9620 575.000 575.0 640
4 0.25 0.45800 0.161000 -10.600 0.9820 482.000 482.0 543
5 0.30 0.37000 0.033700 4.040 0.9830 396.000 395.0 450
6 0.35 0.29100 0.000637 -0.451 0.9830 319.000 319.0 367
7 0.40 0.22300 0.130000 4.460 0.9660 250.000 251.0 286
8 0.45 0.16500 0.150000 4.620 0.9090 191.000 191.0 227
9 0.50 0.11700 0.129000 -3.820 0.9720 141.000 141.0 178
10 0.55 0.07960 0.014000 -0.927 0.9780 100.000 101.0 129
11 0.60 0.05100 0.009720 -0.618 0.9860 68.300 68.0 94
12 0.65 0.03040 0.090100 -1.730 0.9930 43.400 43.0 70
13 0.70 0.01650 0.166000 -1.540 0.9710 25.800 25.5 46
14 0.75 0.00785 0.251000 -1.540 0.9820 13.600 13.0 29
15 0.80 0.00311 0.221000 -3.510 0.2660 6.310 6.0 20
16 0.85 0.00092 0.741000 -2.380 0.9680 2.240 2.0 10
17 0.90 0.00016 0.066100 -2.640 -0.0933 0.478 0.0 4
$cutEstimate
[1] NA
$fitIndices
Cut p.value SFT.R.sq slope. truncated.R.2 mean.k. median.k. max.k.
1 0.10 0.7698749999 0.1290790714 17.9939990 0.97500919 781.016 781.0 825
2 0.15 0.6597824485 0.1133110002 -13.0200677 0.95706793 675.634 676.0 736
3 0.20 0.5554454421 0.1596124481 -11.6012789 0.96210336 575.292 575.0 640
4 0.25 0.4584436539 0.1614028198 -10.5888632 0.98160376 481.764 482.0 543
5 0.30 0.3700831223 0.0337023377 4.0407518 0.98286505 395.628 395.0 450
6 0.35 0.2913484779 0.0006370488 -0.4512638 0.98253470 318.634 319.0 367
7 0.40 0.2228683501 0.1297646029 4.4609108 0.96566340 250.474 251.0 286
8 0.45 0.1648959909 0.1497857583 4.6248899 0.90927965 190.910 191.0 227
9 0.50 0.1173068030 0.1289897777 -3.8171609 0.97150651 140.946 141.0 178
10 0.55 0.0796140219 0.0139719875 -0.9272260 0.97817099 100.418 101.0 129
11 0.60 0.0510032607 0.0097188946 -0.6175709 0.98593371 68.276 68.0 94
12 0.65 0.0303859638 0.0900913765 -1.7328536 0.99324252 43.444 43.0 70
13 0.70 0.0164709788 0.1658053604 -1.5374235 0.97092572 25.808 25.5 46
14 0.75 0.0078523983 0.2505013985 -1.5377157 0.98227575 13.592 13.0 29
15 0.80 0.0031104283 0.2207179089 -3.5116139 0.26596937 6.306 6.0 20
16 0.85 0.0009201220 0.7413987164 -2.3807465 0.96838864 2.242 2.0 10
17 0.90 0.0001599714 0.0660983194 -2.6389376 -0.09327525 0.478 0.0 4
How odd, I have 500 genes and 100 samples so it would be a 100x500 matrix. I've checked on that size and I hit the same issue
Sorry just to add to this, when you increase your example to 100 samples rather than 10 it gives the same error. So it seems its not just down to the number for genes.
> pickHardThreshold(exp_matrix, cutVector = seq(0.1, 0.6, by = 0.1))
Cut p.value SFT.R.sq slope. truncated.R.2 mean.k. median.k. max.k.
1 0.1 0.770 0.032300 5.66 0.1860 6.6 7.0 8
2 0.2 0.555 0.001530 -0.64 0.0390 4.6 5.0 7
3 0.3 0.370 0.000171 -0.14 0.1310 3.0 3.5 4
4 0.4 0.223 0.014800 -1.32 -0.1340 2.2 2.0 4
5 0.5 0.117 0.030900 -1.89 -0.0361 1.0 1.0 3
6 0.6 0.051 0.113000 -3.05 0.7060 0.2 0.0 1
$cutEstimate
[1] NA
$fitIndices
Cut p.value SFT.R.sq slope. truncated.R.2 mean.k. median.k. max.k.
1 0.1 0.76987500 0.0323011476 5.6557635 0.18627071 6.6 7.0 8
2 0.2 0.55544544 0.0015318734 -0.6395621 0.03898946 4.6 5.0 7
3 0.3 0.37008312 0.0001711509 -0.1400856 0.13106798 3.0 3.5 4
4 0.4 0.22286835 0.0148168861 -1.3219051 -0.13430132 2.2 2.0 4
5 0.5 0.11730680 0.0308509930 -1.8857022 -0.03613967 1.0 1.0 3
6 0.6 0.05100326 0.1133166019 -3.0497797 0.70582124 0.2 0.0 1
> pickHardThreshold(exp_matrix, = seq(0.1, 0.7, by = 0.1))
Cut p.value SFT.R.sq slope. truncated.R.2 mean.k. median.k. max.k.
1 0.1 0.7700 0.032300 5.66 0.1860 6.6 7.0 8
2 0.2 0.5550 0.001530 -0.64 0.0390 4.6 5.0 7
3 0.3 0.3700 0.000171 -0.14 0.1310 3.0 3.5 4
4 0.4 0.2230 0.014800 -1.32 -0.1340 2.2 2.0 4
5 0.5 0.1170 0.030900 -1.89 -0.0361 1.0 1.0 3
6 0.6 0.0510 0.113000 -3.05 0.7060 0.2 0.0 1
7 0.7 0.0165 0.113000 -3.05 0.7060 0.2 0.0 1
$cutEstimate
[1] NA
$fitIndices
Cut p.value SFT.R.sq slope. truncated.R.2 mean.k. median.k. max.k.
1 0.1 0.76987500 0.0323011476 5.6557635 0.18627071 6.6 7.0 8
2 0.2 0.55544544 0.0015318734 -0.6395621 0.03898946 4.6 5.0 7
3 0.3 0.37008312 0.0001711509 -0.1400856 0.13106798 3.0 3.5 4
4 0.4 0.22286835 0.0148168861 -1.3219051 -0.13430132 2.2 2.0 4
5 0.5 0.11730680 0.0308509930 -1.8857022 -0.03613967 1.0 1.0 3
6 0.6 0.05100326 0.1133166019 -3.0497797 0.70582124 0.2 0.0 1
7 0.7 0.01647098 0.1133166019 -3.0497797 0.70582124 0.2 0.0 1
If I modify cutVectorI am able to run the function without errors. So i guess the error happen when pickHardThreshold attempt to calculate the scale free topology fit from scaleFreeFitIndex(k = khelp, nBreaks = nBreaks, removeFirst = removeFirst) using a vector k of connectivities which, in your case, seems to be Zero. Funny enough, this does not happen if you use the function pickSoftThreshold.
I would post the same question on bioconductor. WGCNA devs are more active there
Perfect, I will do. Thank you very much for all the help!
Log in to answer this question.
