Molecular term symbols for excited state oxygen
What are the possible term symbols for excited state oxygen with configuration $(1\pi_u)^3(1\pi^*_g)^3$ (all other occupied orbitals are closed shell)? Based on the possible values of $S$ and $\lambda$ it seems that $^1\Delta_u$, $^3\Delta_u$, $^1\Sigma_u$ and $^3\Sigma_u$ are possible.


Term symbols for $(1\pi_u)^3(1\pi_g^*)^3$ - i.e. 1 hole + 1 hole - is equivalent to finding terms for $(1\pi_u)^1(1\pi_g^*)^1$, i.e. 1 electron + 1 electron. I will continue answering the question for the latter case as this is easier to handle.

The fastest way is to use group theory. $O_2$ molecule belongs to $D_\infty^h$ point group. We look at the character table for this group (e.g. here) and perform a direct product of irreducible representations corresponding to our electrons, i.e.:

$\Pi_u \otimes \Pi_g = \Sigma_u^+ \oplus \Sigma_u^- \oplus \Delta_u$

Note this only pertains to spatial symmetry, i.e. orbital angular momentum. It tells us nothing about possible spin states. However, since the two electrons are in different orbitals, all spin combinations are allowed in all cases. The resulting terms are:

TLDR: $^3\Sigma_u^+,$$^3\Sigma_u^-,$$^3\Delta_u,$$^1\Sigma_u^+,$$^1\Sigma_u^-,$$^1\Delta_u$

Your answer is thus correct, the missing piece being the +/- signs of the $\Sigma$ terms. Note that in answering we relied heavily on Pauli exclusion principle not playing a role - since the two electrons occupy different orbitals. Were that not the case, some of the terms would disappear (terms of e.g. $(1\pi_u)^1(2\pi_u)^1$ are $^3\Sigma_g^+,$$^3\Sigma_g^-,$$^3\Delta_g,$$^1\Sigma_g^+,$$^1\Sigma_g^-,$$^1\Delta_g$, while terms of $(1\pi_u)^2$ are a non-trivial subset, specifically $^3\Sigma_g^-,$$^1\Sigma_g^+,$$^1\Delta_g$).


  However, I am not sure in this case how to determine the +/- labels for the $\Sigma$ term symbols.


The following derivations draw heavily from this pdf, read for further details.

We refer to definition of +/- signs as describing the symmetry with respect to reflection in plane that contains both atoms, let's use an arbitrary choice of xz-plane. Wavefunctions symmetric when reflected in this plane are "+" terms, anti-symmetric are "-" terms. However, imagining reflections is not practical. Therefore, we introduce operator notation. The "+" term corresponds to a wave function with eigenvalue of +1 (with respect to $\hat{\sigma}_{xz}$ operator); "-" term corresponds to eigenvalue of -1.

We note:

$$\hat{\sigma}_{xz} \pi_x = \pi_x$$
$$\hat{\sigma}_{xz} \pi_y = - \pi_y$$
$$\pi_1 = -N(\pi_x + i\pi_y)$$
$$\pi_{-1} = N(\pi_x - i\pi_y)$$

, where N is some boring normalization constant. Following up:
$$\hat{\sigma}_{xz} \pi_1 = - \pi_{-1}$$
$$\hat{\sigma}_{xz} \pi_{-1} = - \pi_1$$

We need combinations of configurations (Slater determinants = SD) that are eigenfunctions of $\hat{\sigma}_{xz}$.
E.g.:

$$^1\Sigma^+_u = N(|SD_1\rangle - |SD_2\rangle + |SD_3\rangle - |SD_4\rangle) =$$$$= N(|\pi_{u,1},\bar{\pi}_{g,-1}\rangle - |\bar{\pi}_{u,1}{\pi}_{g,-1}\rangle + |\pi_{u,-1},\bar{\pi}_{g,1}\rangle - |\bar{\pi}_{u,-1}{\pi}_{g,1}\rangle) $$
Let's apply the $\hat{\sigma}_{xz}$ operator:
$$\hat{\sigma}_{xz}{ }^1\Sigma^+_u = \hat{\sigma}_{xz} (N(|\pi_{u,1},\bar{\pi}_{g,-1}\rangle - |\bar{\pi}_{u,1}{\pi}_{g,-1}\rangle + |\pi_{u,-1},\bar{\pi}_{g,1}\rangle - |\bar{\pi}_{u,-1},{\pi}_{g,1}\rangle)) = $$
$$ = N(|-\pi_{u,-1},-\bar{\pi}_{g,1}\rangle - |-\bar{\pi}_{u,-1}{-\pi}_{g,1}\rangle + |-\pi_{u,1},-\bar{\pi}_{g,-1}\rangle - |-\bar{\pi}_{u,1},{-\pi}_{g,-1}\rangle) =$$
$$ = N(|\pi_{u,-1},\bar{\pi}_{g,1}\rangle - |\bar{\pi}_{u,-1}{\pi}_{g,1}\rangle + |\pi_{u,1},\bar{\pi}_{g,-1}\rangle - |\bar{\pi}_{u,1},{\pi}_{g,-1}\rangle) = $$
$$= N(|SD_3\rangle - |SD_4\rangle + |SD_1\rangle - |SD_2\rangle) = { }^1\Sigma^+_u$$

In going from first to second line we apply the operator, i.e. flip "1" subscript to "-1" and vice versa; in addition we multiply each orbital by "-1". Since we have two orbitals this multiplication disappears: ($-1*-1=1$), which is what we do in going from second to third line. In the last line we see that $\hat{\sigma}_{xz}|SD_1\rangle = |SD_3\rangle$;$\hat{\sigma}_{xz}|SD_2\rangle = |SD_4\rangle$ etc.

We could have used a WF in the form e.g.: $WF = |SD_1\rangle + |SD_3\rangle$. Such a WF is indeed an eigenfunction of $\hat{\sigma}_{xz}$. But it is not an eigenfunction of $\hat{S}^2$.

Similarly, a corresponding "-" term is:
$$^1\Sigma^-_u = N(|SD_1\rangle - |SD_2\rangle - |SD_3\rangle + |SD_4\rangle) =$$
$$\hat{\sigma}_{xz}^1\Sigma^-_u = N(|SD_3\rangle - |SD_4\rangle - |SD_1\rangle + |SD_2\rangle) =$$
$$ = -1 * N(- |SD_3\rangle + |SD_4\rangle + |SD_1\rangle - |SD_2\rangle) = -1 * { }^1\Sigma^-_u$$

The four SDs we use here are all eigenfunctions of $\hat{L}_z$ and $\hat{S}_z$ with $M_l = 0$ (i.e. one $\pi_1$ + one $\pi_{-1}$ orbital) and $M_s = 0$ (i.e. one alpha and one beta orbital). There are two more combinations of these SD and these are part of the $^3\Sigma^-_u$ and $^3\Sigma^+_u$ states.

Lastly, non-$\Sigma$ terms have both a "+" and "-" component that flips to each other upon reflection, e.g.: $\hat{\sigma}_{xz} \Delta^+ = \Delta^-$. But since $\Delta^+$ and $\Delta^-$ are different SDs (different configurations), we don't use the label in these cases - they are reserved for $\Sigma$ terms.

The take-away message here is that "+/-" label relates to some kind of symmetry. Similarly to "g/u" label, it only is relevant in some special molecules ("g/u" in molecules with inversion center; "+/-" in diatomic molecules). We use group theory where possible or a very technical procedure when not.


  For ground state oxygen, triplet states pair with - and singlet states with + to give overall antisymmetric wavefunctions. However, in this case, two different orbitals are involved, so does this rule no longer apply?


There is some confusion here. As mentioned above, for ground state ($(1\pi_u)^2$) the terms are indeed $^3\Sigma_g^-$ and $^1\Sigma_g^+$. But there is also the state $^1\Delta_g$ that has no "+/-" label, yet it is still anti-symmetric. The "+/-" label has nothing to do with anti-symmetry. Anti-symmetry is ensured by use of Slater determinants. The "rule" you mentioned is a result applicable to $\pi^2$ (and other doubly degenerate 2 electron) configurations, but is not related to anti-symmetry and gives no predictions for other configurations.