We present the first detailed analysis of the connection between galaxies and their dark matter halos for the unWISE galaxy catalog -- a full-sky, infrared-selected sample built from WISE data, containing over 500 million galaxies. Using unWISE galaxy-galaxy auto-correlation and Planck CMB lensing-galaxy cross-correlation measurements down to 10 arcmin angular scales, we constrain the halo occupation distribution (HOD), a model describing how central and satellite galaxies are distributed within dark matter halos, for three unWISE} galaxy samples at mean redshifts $\bar{z} \approx 0.6$, $1.1$, and $1.5$. We constrain the characteristic minimum halo mass to host a central galaxy, $M_\mathrm{min}^\mathrm{HOD} = 1.83^{+0.41}_{-1.63} \times 10^{12} M_\odot/h$, $5.22^{+0.34}_{-4.80} \times 10^{12} M_\odot/h$, $6.60 ^{+0.30}_{-1.11} \times 10^{13} M_\odot/h$ for the unWISE samples at $\bar{z}\approx 0.6$, $1.1$, and $1.5$, respectively. We find that all three samples are dominated by central galaxies, rather than satellites. Using our constrained HOD models, we infer the effective linear galaxy bias for each unWISE sample, and find that it does not evolve as steeply with redshift as found in previous perturbation-theory-based analyses of these galaxies. We discuss possible sources of systematic uncertainty in our results, the most significant of which is the uncertainty on the galaxy redshift distribution. Our HOD constraints provide a detailed, quantitative understanding of how the unWISE galaxies populate the underlying dark matter halo distribution. These constraints will have a direct impact on future studies employing the unWISE galaxies as a cosmological and astrophysical probe, including measurements of ionized gas thermodynamics and dark matter profiles via Sunyaev-Zel'dovich and lensing cross-correlations.