This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Error in ComBat function: Error in dat[, batch == batch_level] : (subscript) logical subscript too long

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

1 answer

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)

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"

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

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

Log in to answer this question.