Let me first establish what CBoW and skip-gram are supposed to do. You can skip to the next section if you think this is unnecessary.
Background
My understanding is that Word2Vec is a suite of 2 algorithms, continuous bag-of-words (CBoW) and skip-gram neural-network language model (SGNNLM or simply skip-gram), which are both two-layer neural networks. That is, given a vocabulary $V$ and embedding size $H$, they take vectors of size $|V|$ and pass them through one hidden layer and one output layer, with weight matrices $W_1\in \mathbb{R}^{H\times |V|}$ and $W_2\in \mathbb{R}^{|V|\times H}$, much like a small auto-encoder. The output layer is always cited as having softmax activation, and according to this Coursera course, the hidden layer has ReLU activation. I've seen it cited as having "no" activation ($\varphi(x) = x$), though.
The difference between CBoW and skip-gram is supposedly that CBoW "predicts the current word based on the context", and the skip-gram "predicts surrounding words given the current word". That's a literal quote from the original Word2Vec paper, Mikolov 2013. Here's their figure, which looks intuitive, but in reality is notoriously confusing:
Since both are used to move from a one-hot encoding to a smaller embedding (the weights in one or both of the matrices), their input and output format is still one-hot. I know that per training example, CBoW starts out with $C$ one-hot vectors (one for each context word) and just averages them to get a single $|V|$-sized vector. (This is why it's called "bag of words": because averaging is commutative, you lose the order of the context words.) Apart from this tricky many-to-one transformation, the rest of the network makes sense.
Skip-gram has no issue encoding the input, but it has a one-to-many transformation at the end. This is a problem.
The issue
What does skip-gram do in its final layer? I have been scouring various StackExchange sites to find a consensus on this question, and there seem to be fourdistinct camps.
Skip-gram predicts $C$ softmaxes of size $|V|$ using $C$ different matrices $W_{2}^{(1)} \dots W_{2}^{(C)}$. Examples: this article (at least the figures), this article stating "each with its own weight matrix", and this question.
Skip-gram predicts $C$ softmaxes of size $|V|$ using the same matrix $W_2$. Examples: this question and this answer.
Skip-gram predicts $1$ softmax of size $|V|$ using one matrix $W_2$, and the dataset consist of (center word, context word) pairs instead of (center word, all context words). Examples: this tutorial.
Skip-gram is a binary classifier predicting whether a given word is "in context" vs. "out of context" for another word. The dataset consists of pairs of words that either appeared close to each other (positive) or were randomly paired (negative). Examples: this article, this article, and this answer.
Skip-gram is usually explained as having one output matrix, which is hypothesis (2) or (3), yet producing $C$ different words in the output layer, from the same hidden representation, which is hypothesis (1). See this article and the figure in this question. That's a magic trick. Pure fiction.
I find (3) the most reasonable, since it only has two weight matrices and doesn't produce useless copies. In that case, skip-gram actually doesn't predict any words, but rather produces a single $|V|$-sized vector which represents the average one-hot vector of the context words. This single vector can then account for one term in the loss per context word.
Yet, there are some major issues with (3):
You can never reach $0$ loss. Let's say the loss function is log likelihood, and compute the loss for one context window: $$ \ell(\text{word}, \text{context}) = \sum_{i=1}^C \ln P(\text{context}_i \mid \text{word}) $$ Clearly, that probability cannot be $1$ for each context word. In the perfect case, it is $1/C$.
It cannot handle negative examples like (4). Let $y$ be a binary variable that is $1$ when a pair of words is positive and $0$ when it is negative. What is $P(y \mid w_1,w_2)$ given the above model? It can't just be $P(w_2 \mid w_1)$. Why? Look at this binary loss: $$\begin{aligned} \ell(w_1, w_2, y) &= y\ln P(y = 1 \mid w_1, w_2) + (1-y)\ln (1 - P(y = 1 \mid w_1,w_2)) \\ &= y\ln P(w_2 \mid w_1) + (1-y)\ln (1 - P(w_2 \mid w_1)) \end{aligned}$$ As we saw before, best-case scenario, $P(w_2 \mid w_1) = 1/C$ for a positive example, which means $P(y = 0 \mid w_1,w_2) = 1 - 1/C$, clearly way higher. For example: with a context of $C = 4$, the model predicts that a word is context by giving it 25% probability, and accidentally assigns the remaining 75% to all other words, so it always predicts that a word is out of context. You simply cannot combine softmax likelihood with negative sampling. This answer specifically says that softmax likelihood isn't used at all, and instead we use "noise-contrastive estimation".
So, once and for all: does Word2Vec's skip-gram NNLM even produce context words?