When a wavefunction collapses on its own what basis does it project onto? (In Quantum Thermo, why do we privilege energy eigenstates?) Obviously when an experimenter measures it, it projects onto the basis that the experimenter was trying to measure. But this process occurs on it's own, all the time. 
Obviously the answer has to do with the exact situation in which it occurs, but what I'm really struggling with is conceptualizing Quantum Thermodynamics because of it. Why do we get our entropy and everything else by counting energy eigenstates? Couldn't the particles in question be in arbitrary combinations of energy eigenstates? If they're entangling and collapsing constantly, why onto the energy eigenbasis?
 A: We actually force the use of which eigenstates by the macroscopic quantities we want to measure. The starting point is the unknown density operator $D$. We then define its entropy
$$\newcommand{\trace}[1]{\,\mathrm{Tr}\left(#1\right)}S = -k \trace{D\log D}.$$
We wish to minimise it under the constraint that the energy is known,
$$\trace{DH} = U,$$
where $H$ is the hamiltonian. This last equation is the way to encode a measure in the density operator formalism. 
Let's continue to the end by writing the density operator that solve the minimisation problem I have just described. This is the Boltzmann-Gibbs
$$D=\frac{1}{Z}\exp\left(-\frac{H}{kT}\right).$$
The key quantity from which pretty much any thermodynamical quantity can be computed is the partition function function $Z$, whose expression comes from the requirement that $D$ shall be normalised, i.e. $\trace{D}=1$,
$$Z=\trace{\exp\left(-\frac{H}{kT}\right)}.$$ 
Now let's explicit the trace to make it clear,
$$\newcommand{\ket}[1]{|{#1}\rangle}\newcommand{\bra}[1]{\langle{#1}|}Z = \sum_\alpha \bra{\alpha}\exp\left(-\frac{H}{kT}\right)\ket{\alpha},$$
for any orthonormal basis $\{\ket{\alpha}\}$. This is usually tractable only if we choose a basis of eigenstate of $H$,
$$H\ket{\alpha}=E_\alpha\ket{\alpha},$$
because then
$$Z=\sum_\alpha\exp\left(-\frac{E_\alpha}{kT}\right).$$
But note that if we impose another macroscopic constraint than the total energy, then we need more than just the eigenstates of $H$. For example, if we impose a macroscopic magnetic moment, then we will typically need the eigenstates for the spins of our system. The story is the same as above: we will have to take the trace of a Boltzmann-Gibbs distribution featuring an exponential of some spin projector instead of the Hamiltonian, and that will be tractable only if we work with the spin eigenstates.
A: It doesn't actually matter. For example, the microcanonical ensemble would have density matrix
$$\rho = \sum_n | E_n \rangle \langle E_n|$$
but this is just the identity, so it is also equal to
$$\rho = \sum_n | v_n \rangle \langle v_n |$$
for any basis $|v_n \rangle$. You might ask, is it really physically equal? Well, the full quantum states might be different, but if the density measurements are the same, then the two are completely indistinguishable as long as you're only measuring the system and not its surroundings too. And this is exactly the point of stat mech / ensembles in general -- even in classical stat mech we never care or talk about how the system gets into equilibrium, just the results of measurements on it.
As another example, for spin systems sometimes people say there's a 50% chance of spin up and a 50% chance of spin down, but this is completely equivalent to a 50% chance of spin left and a 50% chance of spin right. For example, the chance of measuring spin down is 50% in the first example, and (2)(25%) = 50% in the second.
