Selected distributional properties of the maximum likelihood estimator and its z-transformation of three familial correlations (parental, parent-offspring, filial) were investigated numerically for the case of nuclear families with variable sibship size. This investigation was based on six different sets of the three correlations, and four different sample sizes, defining 24 sampling conditions, which were replicated 1,000 times each. It was found that the distributional properties of the correlation estimator are affected by the magnitude of the correlations even in large samples although approximate normality is achieved locally. Fisher's z-transformation, here used only in its interclass form, achieves reduction of skewness, stabilization of variance, and approach to normality already in small samples, except for the filial correlation (where it may be deemed inappropriate) in smaller samples. For both the correlation estimator and its z-transformation, the (estimated) relative efficiency was shown to be high (better than 90% in most sampling conditions), suggesting that the estimated minimum variance bound is a satisfactory estimator of the sampling variance. It is concluded that the maximum likelihood estimation of familial correlations under variable sibship size is feasible and, when prudently applied, especially in the form of their z-transformations, provides an appropriate method in analyses of family studies.