On compting the dirichlet inverse $f^{-1}(pqrs)$ where I assumed $(p, q, r, s)$ to be arbitrary primes, I arrived at the formula:
$(p, q, r, s)$
$$ f^{-1}(pqrs) = 24f(p)f(q)f(r)f(s)/f(1)^5 + \\ -6f(p)f(q)f(rs)/f(1)^4 + -6f(p)f(qr)f(s)/f(1)^4 + -6f(pq)f(r)f(s)/f(1)^4 + \\ -6f(pr)f(q)f(s)/f(1)^4 + -6f(p)f(qs)f(r)/f(1)^4 + -6f(ps)f(q)f(r)/f(1)^4 + \\ \mathbf{2f(p)f(qrs)/f(1)^3} + 2f(pq)f(rs)/f(1)^3 + 2f(pqr)f(s)/f(1)^3 + \\ 2f(pqs)f(r)/f(1)^3 + \mathbf{2f(pr)f(qs)/f(1)^3} + 2f(prs)f(q)/f(1)^3 + \\ 2f(ps)f(qr)/f(1)^3 + \\ -1f(pqrs)/f(1)^2 $$
This was done using code so I'm fairly sure it's correct. What I don't understand is what this is counting. Why are there terms that have both groupings of 2-and-2, as well as groupings of 1-and-3? What is the genral pattern here? Clearly $n!$ is involved, as well as what I thought were ${k}\choose{r}$, but that cannot be since we have groups of unequal size.
Can someone please explain to me what is going on, a proof of what this is counting, and what the general formula is for primes $p_0 p_1 \dots p_n$?
Formal details
the Dirichlet convolution of two functions is defined as:
$$ f, g: \mathbb Z \rightarrow \mathbb C \qquad (f \star g)(n) \equiv \sum_{d|n} f(d) g\left(\frac{n}{d}\right) $$
We can define the inverse function of $f$ with respect to the Dirichlet convolution as the unique function $f^{-1}$ that satisfies the identity $$ (f \star f^{-1})(n) \equiv \begin{cases} 1 & \text{if $n=1$} \\ 0 &\text{otherwise} \end{cases} $$
We can compute this recursively as:
\begin{align*} f^{-1}(1) &= \frac{1}{f(1)} \\ f^{-1}(n) &= \frac{-1}{f(1)} \sum_{d|n~d<n}f\left(\frac{n}{d}\right)f^{-1}(d) \end{align*}
This allows one to symbolically compute $f^{-1}(n)$ for choices of primes. I did this for primes called $p, q, r, s$. The results are:
$(p, q)$
$$ f^{-1}(p, q) = 2f(p)f(q)/f(1)^3 + -1f(pq)/f(1)^2 $$
$(p, q, r)$
$$ f(p, q, r) = -6f(p)f(q)f(r)/f(1)^4 + \\ 2f(p)f(qr)/f(1)^3 + 2f(pq)f(r)/f(1)^3 + 2f(pr)f(q)/f(1)^3 + \\ -1f(pqr)/f(1)^2 $$
Upto this point, one can harbour hope that it is indeed the binomial coefficients, but the part posted above shows that this is not the case, hence the question.
What the coefficients count are the number of ordered factorisations of $n$ corresponding to a multiplicative partition of $n$.
An ordered factorisation of $n$ with $k$ factors is a(n ordered) $k$-tuple $(d_1, \dotsc, d_k)$ of integers $d_{\kappa} \geqslant 2$ such that $$n = \prod_{\kappa = 1}^{k} d_{\kappa}\,.$$ Let us denote the set of all ordered factorisations of $n$ with $k$ factors by $OF(n,k)$. It is easy to see that $OF(n,k) = \varnothing$ for $k > \Omega(n)$ (the number of prime factors of $n$ counted with multiplicity), and $OF(n,0) = \varnothing$ for $n > 1$, while $OF(1,0) = \{()\}$ is the singleton set containing the $0$-tuple $()$.
Now let us first find an expression for $f^{-1}(n)$ in terms of ordered factorisations of $n$. Iterating from the fundamental $$f^{-1}(n) = -\frac{1}{f(1)}\sum_{\substack{d \mid n \\ d > 1}} f(d)f^{-1}\biggl(\frac{n}{d}\biggr) \tag{1}$$ for $n > 1$, we obtain $$f^{-1}(n) = \sum_{k = 0}^{\Omega(n)} \frac{(-1)^k}{f(1)^{k+1}} \sum_{(d_1,\dotsc, d_k) \in OF(n,k)} \prod_{\kappa = 1}^{k} f(d_{\kappa})\,. \tag{2}$$ We prove $(2)$ by induction.
For $n = 1$, $\Omega(n) = 0$ and the right hand side of $(2)$ evaluates to $1/f(1)$, which indeed is $f^{-1}(1)$.
Let now $n > 1$, and assume $(2)$ already established for all $m < n$. By $(1)$, we have $$f^{-1}(n) = -\frac{1}{f(1)}\sum_{\substack{d \mid n \\ d > 1}} f(d)f^{-1}\biggl(\frac{n}{d}\biggr)\,.$$ There the arguments of $f^{-1}$ on the right hand side are all $< n$, hence by the induction hypothesis \begin{align} f^{-1}(n) &= -\frac{1}{f(1)}\sum_{\substack{d \mid n \\ d > 1}} f(d) \sum_{k = 0}^{\Omega(n/d)} \frac{(-1)^k}{f(1)^{k+1}}\sum_{OF(n/d,k)} \prod_{\kappa = 1}^{k} f(d_{\kappa}) \\ &= \sum_{\substack{d\mid n \\ d > 1}} \sum_{k = 0}^{\Omega(n/d)} \frac{(-1)^{k+1}}{f(1)^{k+2}}\sum_{OF(n/d,k)} f(d)\prod_{\kappa = 1}^{k} f(d_{\kappa})\,. \end{align} Now, since $OF(n/d,k) = \varnothing$ for $k > \Omega(n/d)$ we can extend the second sum to $\Omega(n)-1 \geqslant \Omega(n/d)$ without changing the value, and then change the order of summation to obtain \begin{align} f^{-1}(n) &= \sum_{\substack{d\mid n \\ d > 1}} \sum_{k = 0}^{\Omega(n)-1} \frac{(-1)^{k+1}}{f(1)^{k+2}}\sum_{OF(n/d,k)} f(d)\prod_{\kappa = 1}^{k} f(d_{\kappa}) \\ &= \sum_{k = 0}^{\Omega(n)-1} \frac{(-1)^{k+1}}{f(1)^{k+2}} \sum_{\substack{d\mid n \\ d > 1}}\sum_{OF(n/d,k)} f(d)\prod_{\kappa = 1}^{k} f(d_{\kappa})\,. \end{align} The map $(d,(d_1,\dotsc,d_k)) \mapsto (d,d_1,\dotsc,d_k)$ gives a bijection from $$\bigcup_{\substack{d \mid n \\ d > 1}} \{d\} \times OF(n/d,k)$$ to $OF(n,k+1)$, thus \begin{align} f^{-1}(n) &= \sum_{k = 0}^{\Omega(n)-1} \frac{(-1)^{k+1}}{f(1)^{k+2}} \sum_{\substack{d\mid n \\ d > 1}}\sum_{OF(n/d,k)} f(d)\prod_{\kappa = 1}^{k} f(d_{\kappa}) \\ &= \sum_{k = 0}^{\Omega(n)-1} \frac{(-1)^{k+1}}{f(1)^{k+2}} \sum_{OF(n,k+1)} \prod_{\kappa = 1}^{k+1}f(d_{\kappa}) \\ &= \sum_{k = 1}^{\Omega(n)} \frac{(-1)^k}{f(1)^{k+1}} \sum_{OF(n,k)} \prod_{\kappa = 1}^{k} f(d_{\kappa}) \\ &= \sum_{k = 0}^{\Omega(n)} \frac{(-1)^k}{f(1)^{k+1}} \sum_{OF(n,k)} \prod_{\kappa = 1}^{k} f(d_{\kappa}) \end{align} where we reindexed in the penultimate line and used $OF(n,0) = \varnothing$ for $n > 1$ in the last. Thus $(2)$ also holds for $n$, and our induction proof is complete.
The product $\prod_{\kappa = 1}^{k} f(d_{\kappa})$ of course depends only on how often each divisor of $n$ occurs in the $k$-tuple $(d_1, \dotsc, d_k)$ and not on the order in which they occur. Thus we can collect all the $k$-tuples that can obtained from each other by a permutation of the components into one group. Such a group corresponds to a multiplicative partition of $n$.
We can view a multiplicative partition of $n$ as a function $$u \colon D'(n) \to \mathbb{N}\cup \{0\}$$ such that $$\prod_{d \in D'(n)} d^{u(d)} = n\,,$$ where $D'(n) = \{ d : d \mid n, d > 1\}$. Let $MP(n)$ be the set of all multiplicative partitions of $n$, and for $u \in MP(n)$ let $$\sigma(u) = \sum_{d \in D'(n)} u(d)\,,$$ that is, $\sigma(u)$ is the number of factors in the partition. The number of ordered factorisations of $n$ corresponding to the multiplicative partition $u$ is $$N(u) := \frac{\sigma(u)!}{\prod_{d \in D'(n)} u(d)!}$$ and therefore we can express the formula for $f^{-1}(n)$ in terms of multiplicative partitions of $n$ as $$f^{-1}(n) = \sum_{u \in MP(n)} \frac{(-1)^{\sigma(u)}}{f(1)^{\sigma(u)+1}}\cdot N(u) \prod_{d \in D'(n)} f(d)^{u(d)}\,. \tag{3}$$ If $n$ is squarefree, a multiplicative partition of $n$ cannot have two equal factors, and thus $N(u)$ is simply $\sigma(u)!$ in that case.