Motivated by the negative binomial distribution that is widely used in ecological and biological applications, we develop a new clustered occurrence model where presence/absence data are modelled under a multivariate negative binomial framework. The proposed model allows the estimation of negative binomial parameters from clustered occurrence data. First, we provide conditions to show the existence of maximum likelihood estimates when cluster sizes are homogeneous and equal to two or three. We also consider a composite likelihood approach which allows for additional robustness and flexibility in fitting a broader class of multivariate Bernoulli models for clustered occurrence data. The proposed method is evaluated in a simulation study and demonstrated using forest plot data from the Center for Tropical Forest Science. Finally, we present several examples using multiple visit occupancy data to illustrate the difference between the proposed model and those of N-mixture models.