Problems with GLMM models to assess nested Count data in R

=======

I want to understand what factors better explain the implementation of a directive. My data set includes Three Fixed variables and one nested. Fixed variables are Biogeographic Region with 2 Levels, Criteria (16 levels) and Groups (11 levels). Nested variable is the Marine Sub-Unit (6 levels and nested within Biog-Reg). Each Marine Sub-unit includes counts concerning Criteria and Group. All my factors are categorical. An example of the data set is the following:

Biog_Reg Marine.Sub.Unit Criteria Group Count

1 Celtic Seas FR_Celtic Seas 1.1 Benthic habitats 1

2 Celtic Seas FR_Celtic Seas 1.1 Fish 15

3 Celtic Seas FR_Celtic Seas 1.1 Marine mammals 6

4 Celtic Seas FR_Celtic Seas 1.1 Marine turtles 6

5 Celtic Seas FR_Celtic Seas 1.1 Pelagic habitats 0

6 Celtic Seas FR_Celtic Seas 1.1 Plankton 0

7 Celtic Seas FR_Celtic Seas 1.1 Rock & Biogenic Reef 0

8 Celtic Seas FR_Celtic Seas 1.1 Seabirds 14

9 Celtic Seas FR_Celtic Seas 1.1 Sedimentary habitat 1

10 Celtic Seas FR_Celtic Seas 1.2 Fish 0

... with 632 more rows

Since I am addressing nested count data and I want to access all proper interactions, I am usign GLMM to perform the analysis. I have tested several models but the most complete includes the following design:

Model3 <- glmer (Count ~ Biog_Reg + Criteria + Group + Biog_Reg: Criteria+ Biog_Reg:Group + (1|Marine.Sub.Unit/Biog_Reg), family = poisson (), data = sum)

(I understand that going instead for the following model would be easier but it simply wouldn't run: Model3 <- glmer (Count ~ Biog_Reg * Criteria * Group + (1|Marine.Sub.Unit/Biog_Reg), family = poisson (), data = sum) )

However, the model takes a very long time to resolve and when it does, provides the following errors:

*"Correlation matrix not shown by default, as p = 54 > 12. Use print(x, correlation=TRUE) or vcov(x) if you need it convergence code: 0 unable to evaluate scaled gradient Model failed to converge: degenerate Hessian with 2 negative eigenvalues failure to converge in 10000 evaluations Warning messages: 1: In vcov.merMod(object, use.hessian = use.hessian) : variance-covariance matrix computed from finite-difference Hessian is not positive definite or contains NA values: falling back to var-cov estimated from RX"

I am relatively new in GLM in general, so I have several questions: 1. Do I have the type of model and model syntaxe correct? My variable is naturally nested, i.e. each Marine Sub-unit only belong to one of the two Biogeogrpahic Region, so should a simpler anova solve this?

2. In spite of the errors, the model provides results for each level of each factor individually. If the model is correct, how can I obtain the results for the factor?

3.Should I nest Criteria and Group within Marine Sub-Unit, evend though they are not ramdom? If so, how is the sintaxe for it?

I have already read a huge amount of information but considering

Schielzeth, H. and Nakagawa, S. (2013), Nested by design: model fitting and interpretation in a mixed model era.

I should go for mixed models and GLMM with Poisson distribution should handle count data. I have encountered many resembling questions but none specifically answered my doubts.

Help!!!