Error in dimnames(x) <- dn : length of 'dimnames' [2] not equal to array extent

[gabrielet]

I am having an issue which I am not able to solve. The error is located in the iter_remle() function and its due to wrong colum names assignment.

I am running ANCOMBC_2.1.1. I can't replicate the error and this is what I have, so far:

> # load raw counts from phyloseq object
> phylo_data <- readRDS(paste0(save_obj, "phylo_data_", exp_name, ".Rds"))
> sample_data(phylo_data)$Treatment <- as.factor(sample_data(phylo_data)$Treatment)
> sample_data(phylo_data)$TimePoint <- as.factor(sample_data(phylo_data)$TimePoint)
> sample_data(phylo_data)$CollectionSite <- as.factor(sample_data(phylo_data)$CollectionSite)
> sample_data(phylo_data)$DNA_RNA <- as.factor(sample_data(phylo_data)$DNA_RNA)

> length(sample_data(phylo_data)$DNA_RNA)
[1] 72
> length(sample_data(phylo_data)$CollectionSite)
[1] 72
> length(sample_data(phylo_data)$TimePoint)
[1] 72
> length(sample_data(phylo_data)$Treatment)
[1] 72
> phylo_data
phyloseq-class experiment-level object
otu_table()   OTU Table:         [ 2206 taxa and 72 samples ]
sample_data() Sample Data:       [ 72 samples by 39 sample variables ]
tax_table()   Taxonomy Table:    [ 2206 taxa by 6 taxonomic ranks ]
> # set seed
> set.seed(131)
> # run ancombc2
> ancom_output_xna <- ancombc2(data=phylo_data,
+                          assay_name="counts",
+                          tax_level="Phylum",
+                          fix_formula="DNA_RNA",
+                          rand_formula="1|CollectionSite")
Error in dimnames(x) <- dn :
  length of 'dimnames' [2] not equal to array extent
In addition: There were 50 or more warnings (use warnings() to see the first 50)

Here you can find the warnings which, I believe, are not really worrisome. They are thrown by the reml fits and the iterative reml steps, within iter_remle():

> warnings()
Warning messages:
1: Model failed to converge with 1 negative eigenvalue: -3.6e+00
2: Model failed to converge with 1 negative eigenvalue: -2.6e-01
3: Model failed to converge with 1 negative eigenvalue: -5.5e+00
4: Model failed to converge with 1 negative eigenvalue: -1.0e+00
5: Model failed to converge with 1 negative eigenvalue: -8.0e-04
6: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv,  ... :
  Model failed to converge with max|grad| = 0.0103023 (tol = 0.002, component 1)
7: Model failed to converge with 1 negative eigenvalue: -3.1e+00
8: Model failed to converge with 1 negative eigenvalue: -4.4e+00
9: Model failed to converge with 1 negative eigenvalue: -9.7e-02
10: Model failed to converge with 1 negative eigenvalue: -4.2e-02
11: Model failed to converge with 1 negative eigenvalue: -6.5e-01
12: Model failed to converge with 1 negative eigenvalue: -7.7e+00
13: Model failed to converge with 1 negative eigenvalue: -2.6e+00
14: Model failed to converge with 1 negative eigenvalue: -4.2e+00
15: Model failed to converge with 1 negative eigenvalue: -2.0e-01
16: Model failed to converge with 1 negative eigenvalue: -4.0e+01
17: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv,  ... :
  unable to evaluate scaled gradient
18: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv,  ... :
  Model failed to converge: degenerate  Hessian with 1 negative eigenvalues
19: Model failed to converge with 1 negative eigenvalue: -9.7e-04
20: Model failed to converge with 1 negative eigenvalue: -7.1e-01
21: Model failed to converge with 1 negative eigenvalue: -4.3e-04
22: Model failed to converge with 1 negative eigenvalue: -2.0e-03
23: Model failed to converge with 1 negative eigenvalue: -8.9e-02
24: Model failed to converge with 1 negative eigenvalue: -4.0e-04
25: Model failed to converge with 1 negative eigenvalue: -6.3e-01
26: Model failed to converge with 1 negative eigenvalue: -6.7e+00
27: Model failed to converge with 1 negative eigenvalue: -1.5e-02
28: Model failed to converge with 1 negative eigenvalue: -3.5e-02
29: Model failed to converge with 1 negative eigenvalue: -8.0e-02
30: Model failed to converge with 1 negative eigenvalue: -3.7e+00
31: Model failed to converge with 1 negative eigenvalue: -6.0e-02
32: Model failed to converge with 1 negative eigenvalue: -4.1e-01
33: Model failed to converge with 1 negative eigenvalue: -3.9e+00
34: Model failed to converge with 1 negative eigenvalue: -4.7e-01
35: Model failed to converge with 1 negative eigenvalue: -4.0e-03
36: Model failed to converge with 1 negative eigenvalue: -1.8e+00
37: Model failed to converge with 1 negative eigenvalue: -2.0e-03
38: Model failed to converge with 1 negative eigenvalue: -9.3e-01
39: Model failed to converge with 1 negative eigenvalue: -1.1e-02
40: Model failed to converge with 1 negative eigenvalue: -2.0e-01
41: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv,  ... :
  unable to evaluate scaled gradient
42: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv,  ... :
  Model failed to converge: degenerate  Hessian with 1 negative eigenvalues
43: Model failed to converge with 1 negative eigenvalue: -4.1e-03
44: Model failed to converge with 1 negative eigenvalue: -5.0e+00
45: Model failed to converge with 1 negative eigenvalue: -1.6e-03
46: Model failed to converge with 1 negative eigenvalue: -4.2e-01
47: Model failed to converge with 1 negative eigenvalue: -7.1e-03
48: Model failed to converge with 1 negative eigenvalue: -2.0e+01
49: Model failed to converge with 1 negative eigenvalue: -7.2e-04
50: Model failed to converge with 1 negative eigenvalue: -3.3e+00

The error, on the other hand, is thrown when trying to assign colnames() to _varhat:

> colnames(var_hat) = fix_eff
Error in dimnames(x) <- dn :
  length of 'dimnames' [2] not equal to array extent

> fix_eff
[1] "(Intercept)" "DNA_RNARNA"
> colnames(var_hat)
[1] "(Intercept)"

For now, I am unable to go further that this...but I am investigating.

[gabrielet]

ok, I think i found the problem which is related to this post.

The problem is that the random formula requires the parenthesis, apparently. If i use the command:

> ancom_output_xna <- ancombc2(data=phylo_data,
+                          assay_name="counts",
+                          tax_level="Phylum",
+                          fix_formula="DNA_RNA",
+                          rand_formula="(1|CollectionSite)",
+                          n_cl = 8)

I don't get any problem, even no warnings.