The number of permutations of $p$ without any run of length $r$ (or greater) is$$\int_0^\infty e^{-t}\cdot\prod_i q_{r,n_i}\!(t)\cdot\mathrm{d}t$$ where $q_{r,n}(t)$ is a polynomial in $t$ defined by $$\sum_{n=0}^\infty q_{r,n}(t)x^n=\exp\left(\frac{t(x-x^r)}{1-x^r}\right)$$

For example, $$ q_{2,1}(t)=t\\ q_{2,2}(t)=\tfrac12t^2-t\\ q_{2,3}(t)=\tfrac16t^3-t^2+t\\ q_{2,4}(t)=\tfrac1{24}t^4-\tfrac12t^3+\tfrac32t^2-t\\ q_{2,5}(t)=\tfrac1{120}t^5-\tfrac16t^4+t^3-2t^2+t$$ $$q_{3,1}(t)=t\\ q_{3,2}(t)=\tfrac12t^2\\ q_{3,3}(t)=\tfrac16t^3-t\\ q_{3,4}(t)=\tfrac1{24}t^4-t^2+t\\ q_{3,5}(t)=\tfrac1{120}t^5-\tfrac12t^3+t^2$$

To compute the integral, expand the product and use the fact that $$\int_0^\infty e^{-t}\cdot t^n\cdot\mathrm{d}t=n!$$

For an explanation of how this works, I can't think of a better introduction than Jair Taylor's paper, "Counting Words with Laguerre Series".

Edit. To make it clear: the only hard part is multiplying together the polynomials. Once they have been expanded fully, the integration is easy. Cross out the integral sign, the $e^{-t}$ and $dt$, and replace every $t^n$ with $n!$ $$\require{cancel}\cancel{\int_0^\infty e^{-t}}\cdot a\cdot t^n\cdot\cancel{\mathrm{d}t}=a\cdot n!$$

Edit 2. The coefficients of $q_{r,n}(t)$ can be derived from the $n$th row of the following table: $$\begin{array}{c|ccccc} &1&t&\tfrac1{2!}t^2&\tfrac1{3!}t^3&\ldots\\ \hline x^0&1&0&0&0&\cdots\\ x^1&0&1&0&0&\ddots\\ x^2&0&&1&0&\ddots\\ x^3&0&&&1&\ddots\\ x^4&\vdots&&&&\ddots \end{array}$$ To continue the table for a given $r$, fill the blank spaces, column by column from left to right, using the rule: $$T(n,c)=\sum_{i=0}T(n-1-ri,c-1)-T(n-r-ri,c-1)$$


Note: Here's an another approach based upon example III.24 Smirnov words from Analytic Combinatorics by Philippe Flajolet and Robert Sedgewick.

In order to answer OPs question we will derive a formula from scratch.

At the end I have added a few notes indicating some common aspects of this answer and the answer from @AndrewWoods regarding the highly instructive paper "Counting Words with Laguerre Series" by the MSE member Jair Taylor which is referenced there.

$$$$

In the following it is convenient to rephrase OPs question in terms of words and letters.

So, instead of $n$ marbles of $c$ colors ($c\leq n$), we consider a $c$-ary Alphabet $\mathcal{A}=\{v_1,\ldots,v_c\}$ containing the $c$ letters $v_j, 1\leq j\leq c$ and words of length $n$.

According to OPs Multiset $P=(1^{n_1}2^{n_2}\ldots c^{n_c})$ with $n_1+n_2+\ldots+n_c=n$ and the number of distinct permutations of $P$ $$|\mathcal{S}_{P}|=\binom{n}{n_1,n_2,\ldots,n_c}=\frac{n!}{n_1!n_2!\cdots n_c!}$$ we are looking for a method to generate words of length $n$ from the multiset $W=(v_1^{n_1}v_2^{n_2}\ldots v_c^{n_c})$ with $n_1+n_2+\ldots+n_c=n$. And we ask:

How many words of the multiset $W$ have a run of length (at least) $k$?


We start with

Some Basics: Generating words

We want to use ordinary generating functions and start with some preliminary considerations.

Let $SEQ(v)=\{\varepsilon,v,vv,vvv,\ldots\}$ denote all words with length $\geq0$ containing only the letter $v$. The empty word is denoted with $\varepsilon$. The corresponding generating function is $$W(v)=(vz)^0+(vz)^1+(vz)^2+\ldots=\frac{1}{1-vz}$$ with the exponent of $z^n$ marking the length $n$ of a word of $v$'s and the coefficient of $z^n$ marking the number of words of length $n$. Let $SEQ^{\leq k}(v)$ be the set of all words of $v$'s with length $\leq k$. The generating function for $SEQ^{\leq k}(v)$ is $$1+(vz)+(vz)^2+\ldots+(vz)^{k}=\frac{1-(vz)^{k+1}}{1-(vz)}$$ Without loss of information and in line with Flajolet we may skip $z$ and we will simply write: $$W^{\leq k}(v)=1+v+v^2+\ldots+v^{k}=\frac{1-v^{k+1}}{1-v}$$.

Now, all words containing letters from the alphabet $\mathcal{A}=\{v_1,\ldots,v_c\}$ can be described as

$$W(v_1,\ldots,v_c)=\frac{1}{1-(v_1+\cdots+v_c)}$$

Smirnov Words: Words without runs

This is an interesting twist. Instead of looking for words having a specified run length, we consider Smirnov Words. These are words having no consecutive equal letters. They can be related to unconstrained words and vice versa in the following way:

We observe: If an unconstrained word is given, we can collaps each run of consecutive letters into a single letter associating a Smirnov word this way. Conversely, starting from a Smirnov word and substituting each letter by a sequence of length $\geq 1$ of this letter we can get all unconstrained words.

So, arbitrary words are derived from Smirnov words by a simultaneous substitution:

$$\mathcal{W}=\mathcal{S}[v_1\mapsto SEQ^{\geq 1}(v_1),\ldots,v_c\mapsto SEQ^{\geq 1}(v_c)]$$

with the generating functions

\begin{align*} W(v_1,\ldots,v_c)=S\left(\frac{v_1}{1-v_1},\ldots,\frac{v_c}{1-v_c}\right).\tag{1} \end{align*}

Observe, that relation (1) determines the generating functions $S$ for Smirnov words implicitely. Now, the inverse function of $\frac{v}{1-v}$ is $\frac{v}{1+v}$ and so we can find the solution:

\begin{align*} S(v_1,\ldots,v_c)=W\left(\frac{v_1}{1+v_1},\ldots,\frac{v_c}{1+v_c}\right)=\left(1-\sum_{j=1}^{c}\frac{v_j}{1+v_j}\right)^{-1}\tag{2} \end{align*}

Words with runs of length $< k$

Now, starting with a Smirnov word we can build from it words with runs of length $<k$ by substituting $$v_j\mapsto v_j+\cdots+v_j^{k-1}$$

Let $\mathcal{S}^{<k}(v_1,\ldots,v_c)$ denote the words containing the letters $v_1,\ldots,v_c$ and with runs of length $<k$. The corresponding generating function is

\begin{align*} S^{<k}(v_1,\ldots,v_c)&=\left(1-\sum_{j=1}^{c}\frac{v_j+\cdots+v_j^{k-1}}{1+v_j+\cdots+v_j^{k-1}}\right)^{-1} =\left(1-\sum_{j=1}^{c}v_j\frac{1-v_j^{k-1}}{1-v_j^k}\right)^{-1} \end{align*}

We introduce $z$ again and denote via the coefficient of operator $[z^n]$ the number of words of length $n$ in $S^{<k}$.

\begin{align*} [z^n]S^{<k}(zv_1,\ldots,zv_c)=[z^n]\left(1-\sum_{j=1}^{c}v_jz\frac{1-(v_jz)^{k-1}}{1-(v_jz)^k}\right)^{-1} \end{align*}

And the number of permutations of the multiset $W=(v_1^{n_1}v_2^{n_2}\ldots v_c^{n_c})$ with $n_1+n_2+\ldots+n_c=n$ and runs of length (at least) $k$ is

\begin{align*} \frac{n!}{n_1!n_2!\cdots n_c!}-[v_1^{n_1}\ldots v_c^{n_c}z^n]\left(1-\sum_{j=1}^{c}v_jz\frac{1-(v_jz)^{k-1}}{1-(v_jz)^k}\right)^{-1}\tag{3} \end{align*}


Example: We apply the formula (3) to OPs example with 4 marbles and 4 colors

$$$$

Words with length $n=4$, runs with length $k<2$

We observe (with support from Wolfram Alpha)

\begin{align*} [z^4]S^{<2}&(v_1z,v_2z,v_3z,v_4z)\\ &=[z^4]\sum_{l\geq 0}\left(\sum_{j=0}^{4}v_jz\frac{1-v_jz}{1-(v_jz)^2}\right)^l\\ &=\sum_{l\geq 0}[z^4]\left(\sum_{j=1}^{4}\left(v_jz-(v_jz)^2+(v_jz)^3-(v_jz)^4\right)\right)^l+\mathcal{O}(z^5)\\ &=-v_1^4-v_2^4-v_3^4-v_4^4&(l=1)\\ &\qquad+2(v_1+v_2+v_3+v_4)(v_1^3+v_2^3+v_3^3+v_4^3)+(v_1^2+v_2^2+v_3^2+v_4^2)^2&(l=2)\\ &\qquad-3(v_1+v_2+v_3+v_4)^2(v_1^2+v_2^2+v_3^2+v_4^2)&(l=3)\\ &\qquad+(v_1+v_2+v_3+v_4)^4&(l=4)\\ &=2v_1^2v_2^2+2v_1^2v_2^3+\ldots+2v_3^2v_4^2\\ &\qquad+6v_1^2v_2v_3+6v_1^2v_2v_4+\ldots+6v_2^2v_3v_4\tag{4}\\ &\qquad+24v_1v_2v_3v_4 \end{align*}

We can focus on the leftmost summand of each line of (4) and take the coefficients $$[v_1^2v_2^2]=2,\qquad[v_1^2v_2v_3]=6\qquad\text{and}\qquad[v_1v_2v_3v_4]=24$$ We observe that these three coefficients correspond to the configurations $(2,2), (2,1,1)$ and $(1,1,1,1)$ of marbles and derive in accordance with OPs example in the column runs of length at least 2:

\begin{align*} [v_1^2v_2^2]&=2\qquad\ \{v_1v_2v_1v_2,v_2v_1v_2v_1\}\\ [v_1^2v_2v_3]&=6\qquad\ \{v_1v_2v_1v_3,v_1v_3v_1v_2,v_1v_2v_3v_1,v_1v_3v_2v_1,v_2v_1v_3v_1,v_3v_1v_2v_1\}\\ [v_1v_2v_3v_4]&=24\qquad\{v_1v_2v_3v_3,\ldots,v_4v_3v_2v_1\} \end{align*}

\begin{array}{lll} (4)&\frac{4!}{4!}-0&=1\\ (3,1)&\frac{4!}{3!1!}-0&=4\\ (2,2)&\frac{4!}{2!2!}-2&=4\\ (2,1,1)&\frac{4!}{2!1!1!}-6&=6\\ (1,1,1,1)&\frac{4!}{1!1!1!1}-24&=0 \end{array}

$$$$

Words with length $n=4$, runs with length $k<3$

We observe

\begin{align*} [z^4]S^{<3}&(v_1z,v_2z,v_3z,v_4z)\\ &=[z^4]\sum_{l\geq 0}\left(\sum_{j=0}^{4}v_jz\frac{1-(v_jz)^2}{1-(v_jz)^3}\right)^l\\ &=\sum_{l\geq 0}[z^4]\left(\sum_{j=1}^{4}\left(v_jz-(v_jz)^3+(v_jz)^4\right)\right)^l+\mathcal{O}(z^5)\\ &=v_1^4+v_2^4+v_3^4+v_4^4&(l=1)\\ &\qquad-2(v_1+v_2+v_3+v_4)(v_1^3+v_2^3+v_3^3+v_4^3)&(l=2)\\ &\qquad+0&(l=3)\\ &\qquad+(v_1+v_2+v_3+v_4)^4&(l=4)\\ &=2v_1^3v_2+2v_1^3v_3+\ldots+2v_3^3v_4\\ &\qquad+6v_1^2v_2^2+6v_1^2v_3^2+\ldots+6v_3^2v_4^2\\ &\qquad+12v_1^2v_2v_3+12v_1^2v_2v_4+\ldots+12v_2^2v_3v_4\tag{5}\\ &\qquad+24v_1v_2v_3v_4 \end{align*}

We can focus on the leftmost summand of each line of (5) and take the coefficients $$[v_1^3v_2]=2,[v_1^2v_2^2]=6,\qquad[v_1^2v_2v_3]=12\qquad\text{and}\qquad[v_1v_2v_3v_4]=24$$ We observe that these four coefficients correspond to the configurations $(3,1),(2,2), (2,1,1)$ and $(1,1,1,1)$ of marbles and derive in accordance with OPs example in the column runs of length at least 3:

\begin{array}{lll} (4)&\frac{4!}{4!}-0&=1\\ (3,1)&\frac{4!}{3!1!}-2&=2\\ (2,2)&\frac{4!}{2!2!}-6&=0\\ (2,1,1)&\frac{4!}{2!1!1!}-12&=0\\ (1,1,1,1)&\frac{4!}{1!1!1!1}-24&=0 \end{array}

$$$$

Words with length $n=4$, runs with length $k<4$

We observe

\begin{align*} [z^4]S^{<4}&(v_1z,v_2z,v_3z,v_4z)\\ &=[z^4]\sum_{l\geq 0}\left(\sum_{j=0}^{4}v_jz\frac{1-(v_jz)^3}{1-(v_jz)^4}\right)^l\\ &=\sum_{l\geq 0}[z^4]\left(\sum_{j=1}^{4}\left(v_jz-(v_jz)^4\right)\right)^l+\mathcal{O}(z^5)\\ &=-v_1^4-v_2^4-v_3^4-v_4^4&(l=1)\\ &\qquad+0&(l=2,3)\\ &\qquad+(v_1+v_2+v_3+v_4)^4&(l=4)\\ &=4v_1^3v_2+4v_1^3v_3+\ldots+4v_3^3v_4\\ &\qquad+6v_1^2v_2^2+6v_1^2v_3^2+\ldots+6v_3^2v_4^2\\ &\qquad+12v_1^2v_2v_3+12v_1^2v_2v_4+\ldots+12v_2^2v_3v_4\tag{6}\\ &\qquad+24v_1v_2v_3v_4 \end{align*}

We can focus on the leftmost summand of each line of (6) and take the coefficients $$[v_1^3v_2]=4,[v_1^2v_2^2]=6,\qquad[v_1^2v_2v_3]=12\qquad\text{and}\qquad[v_1v_2v_3v_4]=24$$ We observe that these four coefficients correspond to the configurations $(3,1),(2,2), (2,1,1)$ and $(1,1,1,1)$ of marbles and derive in accordance with OPs example in the column runs of length at least 4:

\begin{array}{lll} (4)&\frac{4!}{4!}-0&=1\\ (3,1)&\frac{4!}{3!1!}-4&=0\\ (2,2)&\frac{4!}{2!2!}-6&=0\\ (2,1,1)&\frac{4!}{2!1!1!}-12&=0\\ (1,1,1,1)&\frac{4!}{1!1!1!1}-24&=0 \end{array}


Notes:

  • In order to obtain the generating functions for words with length of runs $<k$ we put the focus on Smirnov words. The same approach is done in Jair Taylor paper "Counting Words with Laguerre Series". Interestingly these words are called Carlitz words (after Leonard Carlitz) in his paper.

  • Flajolet followed the treatment of Combinatorial Enumeration by Ian P Goulden and David M. Jackson. We can find this example in their book in section 2.4.16 (The Smirnov Problem).

  • Based upon (2) we found the generating function for words with runs $<k$ by substituting $v_j \mapsto v_j+\cdots+v_j^{k-1}$ and we obtained

\begin{align*} S^{<k}(v_1,\ldots,v_c)=\left(1-\sum_{j=1}^{c}v_j\frac{1-v_j^{k-1}}{1-v_j^k}\right)^{-1}\tag{7} \end{align*}

If we consider the univariate problem we get the number of $c$-ary words with runs of length $<k$ and according to (7)

\begin{align*} \left(1-cz\frac{1-z^{k-1}}{1-z^k}\right)^{-1}=\frac{1-z^{k}}{1-cz+(c-1)z^{k}} \end{align*}

Jair Taylor also presents this formula in his paper at the end of section 3. While Flajolets example is based upon generating functions is Jair Taylors approach based upon Laplace transforms $\mathcal{L}$ operating on Laguerre series. He so obtains the relation

\begin{align*} \Phi\left(\exp\left(tcz\frac{1-z^{k-1}}{1-z^k}\right)\right)=\frac{1-z^{k}}{1-cz+(c-1)z^{k}} \end{align*} with $\Phi\left(e^{tf}\right)=\frac{1}{1-f}$.