Limited attention from the modelling community has been given to ecological approaches which aim to predict geographical patterns of malaria by accounting for the joint effects of different vectors and environmental drivers. A hierarchical multivariate joint spatial Gaussian generalised linear model was developed to provide joint parameters inference and mapping of counts of Anopheles gambiae, An. funestus, An. nili and malaria incidence collected in an area of Cote d'Ivoire. Variable‐selection methods were applied to select important predictors for each mosquito species and malaria incidence. The proposed joint model led to a general reduction of the variance in the estimates compared to independent modelling. There was high variability in the composition of Anopheles mosquito species in the villages with each species suitability only partly overlapping geographically. Abundances of An. gambiae, An. funestus and An. nili were primarily determined by temperature. None of the species were found as a significant predictor for the others. Anopheles gambiae was the predominant species and only An. gambiae female abundance was an important variable (linear predictor) for malaria incidence. However, the geographic correlation analyses show that the rest of Anopheles species are likely playing a role in malaria suitability. Residuals from the models of mosquito abundance and malaria cases are also correlated with each other and overlapping but in geographic patches, meaning that local drivers of vector‐malaria suitability are still present and not represented by the predictors used in the model. Synthesis and applications: Joint modelling improve predictive estimation compared to individual modelling. The accurate predictions highlighted high diversity in the association between malaria and vector species, with most of the area having more than one species suitability correlated with malaria suitability. These zones are unlikely to benefit from species‐specific interventions. Areas with correlated malaria and vector species suitability residuals contain local information, not included in the model, that requires further investigation. This will identify additional communal malaria and vectors factors that need to be considered for optimal malaria control and elimination strategies since these factors are expected to be linked to the local malaria transmission.