We present a multivariate version of the recently proposed extended Poisson-Tweedie regression models to deal with multiple response variables in the context of count data. Similarly to the one response variable case the count nature of the data is taken into account by means of the power dispersion function associated with the Poisson-Tweedie distribution. The mean structure of each response variable is modelled by means of a link function and a linear predictor. The covariance structure for each response variable is defined in terms of a covariance link function combined with a matrix linear predictor involving known matrices, while the joint covariance matrix for the multivariate count response vector is specified using the generalized kronecker product. This specification provides a flexible and efficient multivariate regression methodology for a comprehensive family of count models including multivariate analagous to the Hermite, Neyman Type A, Pólya–Aeppli, negative binomial and Poisson-inverse Gaussian. We discuss extensions of the orthodox multivariate analysis of variance (MANOVA) for count data. Furthermore, we present the computational implementation in R through the package mcglm. Illustrations include a five response variables regression model for health services care and a bivariate longitudinal data in the context of bushmeat in Pico Basilé, Bioko Island, Equatorial Guinea.