Probability distribution models are essential for the reliability-based analysis and design of offshore structures. The objective of the present research was developing a probability distribution model for the local joint flexibility (LJF) factor, fLJF, in two-planar tubular DK-joints commonly found in jacket-type offshore wind turbines and oil/gas production platforms. A total of 324 finite element (FE) analyses were carried out on 81 FE models of DK-joints subjected to four types of out-of-plane bending (OPB) moment loads. Generated FE models were validated using the available experimental data, FE results, and design formulas. Based on the results of parametric FE study, a sample database was prepared for the fLJF values and density histograms were generated for respective samples based on the Freedman-Diaconis rule. Nine theoretical probability density functions (PDFs) were fitted to the developed histograms and the maximum likelihood (ML) method was applied to evaluate the parameters of fitted PDFs. In each case, the Kolmogorov-Smirnov and chi-squared tests were used to evaluate the goodness of fit. Finally, the lognormal model was proposed as the governing probability distribution function for the fLJF. After substituting the values of estimated parameters, four fully defined PDFs were presented for the fLJF in tubular DK-joints subjected to four types of OPB moment loads.