A first principles study of the dynamics of $^6$Li($^{2}$S) + $^6$Li$^{174}$Yb($^2Sigma^+$)$ to ^6$Li$_2(^1Sigma^+$) + $^{174}$Yb($^1$S) reaction is presented at cold and ultracold temperatures. The computations involve determination and analytic fitting of a three-dimensional potential energy surface for the Li$_2$Yb system and quantum dynamics calculations of varying complexities, ranging from exact quantum dynamics within the close-coupling scheme, to statistical quantum treatment, and universal models. It is demonstrated that the two simplified methods yield zero-temperature limiting reaction rate coefficients in reasonable agreement with the full close-coupling calculations. The effect of the three-body term in the interaction potential is explored by comparing quantum dynamics results from a pairwise potential that neglects the three-body term to that derived from the full interaction potential. Inclusion of the three-body term in the close-coupling calculations was found to reduce the limiting rate coefficients by a factor of two. The reaction exoergicity populates vibrational levels as high as $v=19$ of the $^6$Li$_2$ molecule in the limit of zero collision energy. Product vibrational distributions from the close-coupling calculations reveal sensitivity to inclusion of three-body forces in the interaction potential. Overall, the results indicate that a simplified model based on the long-range potential is able to yield reliable values of the total reaction rate coefficient in the ultracold limit but a more rigorous approach based on statistical quantum or quantum close-coupling methods is desirable when product rovibrational distribution is required.