Isoniazid (INH) is a front-line drug used in the treatment of tuberculosis (TB), a disease that remains a major cause of death worldwide. Isoniazid is a prodrug, requiring activation in the mycobacterial cell by the catalase-peroxidase (CP) enzyme. Recent studies have suggested that acetylation of INH by the arylamine-N-acetyltransferase from Mycobacterium tuberculosis (TBNAT) may be a possible cause of inactivation of the drug thus resulting in resistant strains. In this study, computational techniques were applied to investigate the binding of isoniazid to three TBNAT isoforms: wild type, G68R and L125M. Since there is no experimental structure available, molecular dynamics (MD) simulations were initially used for the refinement of TBNAT homology models. Distinct conformations of the models were selected during the production stage of MD simulations for molecular docking experiments with the drug. Finally, each mode of binding was refined by new molecular MD simulations. Essential dynamics (ED) analysis and linear interaction energy calculations (LIE) were used to evaluate the impact of amino acid substitutions on the structural and binding properties of the enzymes. The results suggest that the wild type and the G68R TBNATs have a similar pattern of affinity to INH. On the other hand, the calculated enzyme-INH dissociation constant (KD) was estimated 33 times lower for L125M isoform in comparison with wild type enzyme. This last finding is consistent with the hypothesis that isolated mutations in the tbnat gene can produce M. tuberculosis strains resistant to isoniazid.