We performed a Monte Carlo power analysis to obtain an indication on the size of the genetic effects detectable in our sample. To this end, we simulated 10 000 samples consisting of 3690 families of various configurations reflecting the unbalanced structure of families included in the analyses, i.e., families consisting of singletons, two parents or families comprising sibships sizes 1–6 with 0, 1 or 2 parents. Genotypes in Hardy–Weinberg equilibrium were generated at a locus with a MAF of 0.5 and explaining 1.5 and 1 % variance in the phenotype. The normally distributed phenotype was simulated conditional on the locus and then dichotomized using a cut-off point corresponding to a z-score of 0.85 to mimic the 20 % prevalence of initiation observed in the NTR sample. The correlations between spouses, full siblings and parent-offspring estimated in our sample equaled 0.39, 0.35 and 0.15, respectively. An α = 1 × 10−8 was used to assess the power to detect association. To model association we used a generalized equations estimation (GEE) procedure with an exchangeable working correlation matrix and a sandwich correction to correct the standard errors for misspecification of the background model (Minica et al. 2014).