Error in ComBat function: Error in dat[, batch == batch_level] : (subscript) logical subscript too long
1
0
Entering edit mode
3.4 years ago
Sib ▴ 60

Hello Biostars

I want to perform a differential expression analysis after batch effect removal between two datasets. But after running the ComBat function I get an Error.

I'll be thankful if you say me what is wrong.

Here are the codes I run:

library(affy)
library(sva)

#### Read datasets
dataset1<-ReadAffy(celfile.path = "Data/GSE32323/")
dataset2<-ReadAffy(celfile.path = "Data/GSE8671/")
gse<-merge(dataset1,dataset2)

#### Background correction, Normalization, Log transform, and extraction of series matrix
rma<-rma(gse)

#### Batch effect removal
meta <- read.delim("Data/Metadata.txt", row.names=1)
design <- model.matrix(~Tissue, meta)
bc <- ComBat(rma, meta$Batch, design)

And the the error is:

Error in dat[, batch == batch_level] : 
  (subscript) logical subscript too long

Also here is the variable meta:

  Tissue     Batch
GSM800742 Normal     1
GSM800743 Normal     1
GSM800744 Normal     1
GSM800745 Normal     1
GSM800746 Normal     1
GSM800747 Normal     1
GSM800748 Normal     1
GSM800749 Normal     1
GSM800750 Normal     1
GSM800751 Normal     1
GSM800752 Normal     1
GSM800753 Normal     1
GSM800754 Normal     1
GSM800755 Normal     1
GSM800756 Normal     1
GSM800757 Normal     1
GSM800758 Normal     1
GSM800759 Cancer     1
GSM800760 Cancer     1
GSM800761 Cancer     1
GSM800762 Cancer     1
GSM800763 Cancer     1
GSM800764 Cancer     1
GSM800765 Cancer     1
GSM800766 Cancer     1
GSM800767 Cancer     1
GSM800768 Cancer     1
GSM800769 Cancer     1
GSM800770 Cancer     1
GSM800771 Cancer     1
GSM800772 Cancer     1
GSM800773 Cancer     1
GSM800774 Cancer     1
GSM800775 Cancer     1
GSM215051 Normal     2
GSM215052 Normal     2
GSM215053 Normal     2
GSM215054 Normal     2
GSM215055 Normal     2
GSM215056 Normal     2
GSM215057 Normal     2
GSM215058 Normal     2
GSM215059 Normal     2
GSM215060 Normal     2
GSM215061 Normal     2
GSM215062 Normal     2
GSM215063 Normal     2
GSM215064 Normal     2
GSM215065 Normal     2
GSM215066 Normal     2
GSM215067 Normal     2
GSM215068 Normal     2
GSM215069 Normal     2
GSM215070 Normal     2
GSM215071 Normal     2
GSM215072 Normal     2
GSM215073 Normal     2
GSM215074 Normal     2
GSM215075 Normal     2
GSM215076 Normal     2
GSM215077 Normal     2
GSM215078 Normal     2
GSM215079 Normal     2
GSM215080 Normal     2
GSM215081 Normal     2
GSM215082 Normal     2
GSM215083 Cancer     2
GSM215084 Cancer     2
GSM215085 Cancer     2
GSM215086 Cancer     2
GSM215087 Cancer     2
GSM215088 Cancer     2
GSM215089 Cancer     2
GSM215090 Cancer     2
GSM215091 Cancer     2
GSM215092 Cancer     2
GSM215093 Cancer     2
GSM215094 Cancer     2
GSM215095 Cancer     2
GSM215096 Cancer     2
GSM215097 Cancer     2
GSM215098 Cancer     2
GSM215099 Cancer     2
GSM215100 Cancer     2
GSM215101 Cancer     2
GSM215102 Cancer     2
GSM215103 Cancer     2
GSM215104 Cancer     2
GSM215105 Cancer     2
GSM215106 Cancer     2
GSM215107 Cancer     2
GSM215108 Cancer     2
GSM215109 Cancer     2
GSM215110 Cancer     2
GSM215111 Cancer     2
GSM215112 Cancer     2
GSM215113 Cancer     2
GSM215114 Cancer     2

And here is the design:

         (Intercept) TissueNormal
GSM800742           1            1
GSM800743           1            1
GSM800744           1            1
GSM800745           1            1
GSM800746           1            1
GSM800747           1            1
GSM800748           1            1
GSM800749           1            1
GSM800750           1            1
GSM800751           1            1
GSM800752           1            1
GSM800753           1            1
GSM800754           1            1
GSM800755           1            1
GSM800756           1            1
GSM800757           1            1
GSM800758           1            1
GSM800759           1            0
GSM800760           1            0
GSM800761           1            0
GSM800762           1            0
GSM800763           1            0
GSM800764           1            0
GSM800765           1            0
GSM800766           1            0
GSM800767           1            0
GSM800768           1            0
GSM800769           1            0
GSM800770           1            0
GSM800771           1            0
GSM800772           1            0
GSM800773           1            0
GSM800774           1            0
GSM800775           1            0
GSM215051           1            1
GSM215052           1            1
GSM215053           1            1
GSM215054           1            1
GSM215055           1            1
GSM215056           1            1
GSM215057           1            1
GSM215058           1            1
GSM215059           1            1
GSM215060           1            1
GSM215061           1            1
GSM215062           1            1
GSM215063           1            1
GSM215064           1            1
GSM215065           1            1
GSM215066           1            1
GSM215067           1            1
GSM215068           1            1
GSM215069           1            1
GSM215070           1            1
GSM215071           1            1
GSM215072           1            1
GSM215073           1            1
GSM215074           1            1
GSM215075           1            1
GSM215076           1            1
GSM215077           1            1
GSM215078           1            1
GSM215079           1            1
GSM215080           1            1
GSM215081           1            1
GSM215082           1            1
GSM215083           1            0
GSM215084           1            0
GSM215085           1            0
GSM215086           1            0
GSM215087           1            0
GSM215088           1            0
GSM215089           1            0
GSM215090           1            0
GSM215091           1            0
GSM215092           1            0
GSM215093           1            0
GSM215094           1            0
GSM215095           1            0
GSM215096           1            0
GSM215097           1            0
GSM215098           1            0
GSM215099           1            0
GSM215100           1            0
GSM215101           1            0
GSM215102           1            0
GSM215103           1            0
GSM215104           1            0
GSM215105           1            0
GSM215106           1            0
GSM215107           1            0
GSM215108           1            0
GSM215109           1            0
GSM215110           1            0
GSM215111           1            0
GSM215112           1            0
GSM215113           1            0
GSM215114           1            0
attr(,"assign")
[1] 0 1
attr(,"contrasts")
attr(,"contrasts")$Tissue
[1] "contr.treatment"
Batch-effect microarray R • 4.4k views
ADD COMMENT
1
Entering edit mode
3.4 years ago

Are you sure that this command is doing what you think:

gse <- merge(dataset1,dataset2)

Can you show the output of

str(rma)
str(meta)
str(design)
ADD COMMENT
1
Entering edit mode

Thanks, I found the problem . I had removed two samples in datasets while had not removed them in metadata.

ADD REPLY
0
Entering edit mode

Thanks to response.

>str(rma)
    'data.frame':   54675 obs. of  96 variables:
     $ GSM800742_chip_array_C06N.H.CEL: num  10.46 5.03 4.57 8.08 3.3 ...
     $ GSM800743_chip_array_C11N.H.CEL: num  10.54 5.64 4.65 7.89 3.45 ...
     $ GSM800744_chip_array_C24N.H.CEL: num  10.41 6.26 5.01 8.31 3.21 ...
     $ GSM800745_chip_array_C27N.H.CEL: num  10.64 6.05 4.72 8.18 3.16 ...
     $ GSM800746_chip_array_C28N.H.CEL: num  10.26 6.49 4.39 8.32 3.14 ...
     $ GSM800747_chip_array_C30N.H.CEL: num  10.35 6.07 5.11 8.57 3.26 ...
     $ GSM800748_chip_array_C31N.H.CEL: num  10.18 6.23 4.88 8.56 3.26 ...
     $ GSM800749_chip_array_C32N.H.CEL: num  10.19 6.63 5 7.5 3.18 ...
     $ GSM800750_chip_array_C33N.H.CEL: num  10.16 5.74 4.94 8.18 3.1 ...
     $ GSM800751_chip_array_C35N.H.CEL: num  10.31 6.53 5 8.25 3.36 ...
     $ GSM800752_chip_array_C36N.H.CEL: num  8.27 6.37 4.95 7.66 3.44 ...
     $ GSM800753_chip_array_C38N.H.CEL: num  10.26 5.92 4.86 7.68 3.2 ...
     $ GSM800754_chip_array_C41N.H.CEL: num  10.27 5.7 4.65 8.59 3.19 ...
     $ GSM800755_chip_array_C42N.H.CEL: num  10.31 5.86 4.95 8.36 3.29 ...
     $ GSM800756_chip_array_C44N.H.CEL: num  10.12 5.95 4.87 7.92 3.2 ...
     $ GSM800757_chip_array_C45N.H.CEL: num  9.87 6.72 4.65 8.52 3.36 ...
     $ GSM800758_chip_array_C47N.H.CEL: num  10.24 6.4 6.24 7.66 3.27 ...
     $ GSM800759_chip_array_C06T.H.CEL: num  9.8 6.87 4.72 7.97 3.56 ...
     $ GSM800760_chip_array_C11T.H.CEL: num  10.2 6.96 4.72 7.71 3.3 ...
     $ GSM800761_chip_array_C24T.H.CEL: num  10.1 6.27 5.33 7.67 3.71 ...
     $ GSM800762_chip_array_C27T.H.CEL: num  10.15 6.29 4.78 8.19 3.16 ...
     $ GSM800763_chip_array_C28T.H.CEL: num  10.27 7 4.76 8.3 3.23 ...
     $ GSM800764_chip_array_C30T.H.CEL: num  9.76 6.91 4.78 7.44 3.2 ...
     $ GSM800765_chip_array_C31T.H.CEL: num  10.05 7.21 4.66 7.89 3.19 ...
     $ GSM800766_chip_array_C32T.H.CEL: num  9.91 7.17 5.94 7.58 3.13 ...
     $ GSM800767_chip_array_C33T.H.CEL: num  10.16 7.1 4.85 7.83 3.31 ...
     $ GSM800768_chip_array_C35T.H.CEL: num  10.47 7.4 4.93 8.4 3.37 ...
     $ GSM800769_chip_array_C36T.H.CEL: num  10.22 7.15 4.77 8.12 3.3 ...
     $ GSM800770_chip_array_C38T.H.CEL: num  9.87 7.05 5.03 7.61 3.26 ...
     $ GSM800771_chip_array_C41T.H.CEL: num  10.08 6.57 5.5 8.25 3.47 ...
     $ GSM800772_chip_array_C42T.H.CEL: num  9.67 7.6 5.04 8.12 3.35 ...
     $ GSM800773_chip_array_C44T.H.CEL: num  9.71 7.48 4.74 8.05 3.31 ...
     $ GSM800774_chip_array_C45T.H.CEL: num  10.26 7.08 4.78 7.91 3.26 ...
     $ GSM800775_chip_array_C47T.H.CEL: num  9.35 5.9 5.33 7.52 3.45 ...
     $ GSM215051.CEL                  : num  10.58 6.82 5.51 8.22 3.05 ...
     $ GSM215052.CEL                  : num  10.57 6.96 4.8 7.92 3.31 ...
     $ GSM215053.CEL                  : num  10.47 6.89 5.04 7.9 3.25 ...
     $ GSM215054.CEL                  : num  10.45 6.38 5.01 7.96 3.23 ...
     $ GSM215055.CEL                  : num  10.35 7.12 4.87 7.87 3.12 ...
     $ GSM215056.CEL                  : num  10.4 6.64 5.03 7.94 3.2 ...
     $ GSM215057.CEL                  : num  10.37 7.08 5.14 7.94 3.09 ...
     $ GSM215058.CEL                  : num  10.35 7.14 4.98 8.02 3.08 ...
     $ GSM215059.CEL                  : num  10.27 6.99 5.07 8.04 3.09 ...
     $ GSM215060.CEL                  : num  10.47 6.73 4.99 8.01 3.13 ...
     $ GSM215061.CEL                  : num  10.05 6.86 4.73 7.73 3.05 ...
     $ GSM215062.CEL                  : num  10.37 7.3 4.98 7.44 3.01 ...
     $ GSM215064.CEL                  : num  10.27 7.28 5.1 8.13 3.04 ...
     $ GSM215065.CEL                  : num  10.36 6.88 5.02 7.85 3.15 ...
     $ GSM215066.CEL                  : num  10.46 6.85 4.85 8.08 3.17 ...
     $ GSM215067.CEL                  : num  10.34 7.27 5.23 8.2 3.07 ...
     $ GSM215068.CEL                  : num  10.61 6.96 5 8.39 3.2 ...
     $ GSM215069.CEL                  : num  10.5 6.88 4.89 7.92 3.04 ...
     $ GSM215070.CEL                  : num  10.51 6.61 5 8.05 3.24 ...
     $ GSM215071.CEL                  : num  10.36 6.63 4.86 7.84 3.26 ...
     $ GSM215072.CEL                  : num  10.18 7.75 4.78 8.08 3.14 ...
     $ GSM215073.CEL                  : num  10.4 6.49 5.05 7.87 3.1 ...
     $ GSM215074.CEL                  : num  10.42 6.83 5.15 8.06 3.15 ...
     $ GSM215075.CEL                  : num  10.45 7.05 4.82 8.04 3.15 ...
     $ GSM215076.CEL                  : num  10.23 6.9 4.95 8.43 3.09 ...
     $ GSM215077.CEL                  : num  10.24 6.72 4.92 7.94 3.31 ...
     $ GSM215078.CEL                  : num  10.44 6.33 4.88 7.94 3.13 ...
     $ GSM215079.CEL                  : num  10.55 6.86 5.03 7.79 3.15 ...
     $ GSM215080.CEL                  : num  10.5 7.02 4.99 7.84 3.13 ...
     $ GSM215081.CEL                  : num  10.34 7 4.93 7.78 3.14 ...
     $ GSM215082.CEL                  : num  10.72 6.71 4.98 8.25 3.21 ...
     $ GSM215083.CEL                  : num  10.66 7.98 5.03 8.1 3.16 ...
     $ GSM215084.CEL                  : num  10.52 7.42 4.64 7.58 3.1 ...
     $ GSM215085.CEL                  : num  10.28 8.05 5.04 7.33 3.1 ...
     $ GSM215086.CEL                  : num  10.38 8.4 5.1 7.84 3.24 ...
     $ GSM215087.CEL                  : num  10.16 8.14 4.52 7.34 3.04 ...
     $ GSM215088.CEL                  : num  10.28 7.18 5.15 8.11 3.33 ...
     $ GSM215089.CEL                  : num  10.62 7.6 4.84 8.21 3.12 ...
     $ GSM215090.CEL                  : num  10.25 8.09 4.88 8.27 3.42 ...
     $ GSM215091.CEL                  : num  10.17 7.75 4.78 8.33 3.26 ...
     $ GSM215092.CEL                  : num  10.27 8.5 4.76 7.95 3.23 ...
     $ GSM215094.CEL                  : num  10.1 8.33 4.82 7.46 3.13 ...
     $ GSM215095.CEL                  : num  10.38 7.79 4.77 7.52 3.1 ...
     $ GSM215096.CEL                  : num  10.38 8.65 4.84 7.94 3.19 ...
     $ GSM215097.CEL                  : num  10.15 8.02 4.95 8.13 3.28 ...
     $ GSM215098.CEL                  : num  10.02 8 4.77 7.93 3.15 ...
     $ GSM215099.CEL                  : num  10.28 8.72 4.89 7.8 3.15 ...
     $ GSM215100.CEL                  : num  10.21 7.75 4.65 8.42 3.28 ...
     $ GSM215101.CEL                  : num  10.5 7.68 5.03 8.26 3.32 ...
     $ GSM215102.CEL                  : num  10.72 7.53 4.99 8.37 3.13 ...
     $ GSM215103.CEL                  : num  10.21 7.54 4.55 7.85 3.17 ...
     $ GSM215104.CEL                  : num  10.36 8.4 4.64 8.16 3.13 ...
     $ GSM215105.CEL                  : num  10.26 7.43 4.98 7.44 2.96 ...
     $ GSM215106.CEL                  : num  10.18 8.35 4.79 7.76 2.88 ...
     $ GSM215107.CEL                  : num  9.87 8.61 4.56 7.7 3.14 ...
     $ GSM215108.CEL                  : num  10.37 7.57 4.79 8.1 3.35 ...
     $ GSM215109.CEL                  : num  10.53 7.48 4.94 7.7 3.38 ...
     $ GSM215110.CEL                  : num  10.3 8.04 4.71 8.09 3.2 ...
     $ GSM215111.CEL                  : num  10.45 7.84 5.06 7.75 3.13 ...
     $ GSM215112.CEL                  : num  10.35 8.19 4.93 7.57 3.05 ...
     $ GSM215113.CEL                  : num  10.37 8.97 4.79 7.42 3.11 ...
     $ GSM215114.CEL                  : num  10.52 7.59 4.87 8.03 3.06 ...

    > str(meta)
    'data.frame':   98 obs. of  2 variables:
     $ Tissue: chr  "Normal" "Normal" "Normal" "Normal" ...
     $ Batch : int  1 1 1 1 1 1 1 1 1 1 ...


    > str(design)
     num [1:98, 1:2] 1 1 1 1 1 1 1 1 1 1 ...
     - attr(*, "dimnames")=List of 2
      ..$ : chr [1:98] "GSM800742" "GSM800743" "GSM800744" "GSM800745" ...
      ..$ : chr [1:2] "(Intercept)" "TissueNormal"
     - attr(*, "assign")= int [1:2] 0 1
     - attr(*, "contrasts")=List of 1
      ..$ Tissue: chr "contr.treatment"
ADD REPLY
1
Entering edit mode

Hi again, therein lies the problem. rma has 96 columns, while meta and design have 98 rows. These need to match in number and order.

I spotted 2 discrepancies already, i.e., GSM215063 and GSM215093. These are missing from rma

ADD REPLY

Login before adding your answer.

Traffic: 2779 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6