Applying a mask to fmri data

I loaded an fMRI data as 4D data and I want to mask the gray matter. The gray matter mask is probabilistic and it has different dimensions than the fMRI data. So I did the following:

    mask = nib.load(GM_mask)
    mask = mask.get_fdata()
    metadata_holder = nib.load(bold_fname)
    masking_threshold = .5
    mask[mask >= masking_threshold] = 1
    mask[mask < masking_threshold] = 0
    # resample mask to have the same dim. of the fMRI
    mask = zoom(mask, (metadata_holder.shape[0] / mask.shape[0], metadata_holder.shape[1] / mask.shape[1], metadata_holder.shape[2] / mask.shape[2]))
    data = metadata_holder.get_fdata()

    data = mask[:, :, :, np.newaxis] * data

The last line will just element-wise multiply the resampled mask with data. However, I am concerned about the affine transformation, where the voxels spatial locations in the mask and gray matters are different. So, is my approach correct? Or should I stick with the API that nibabel offers?

data = nib.masking.apply_mask(metadata_holder, GM_mask)