首先扔一个结论：记最终答案是 ansans，另外 S(n,m)S(n,m) 表示前 nn 个 自然数 的 mm 次幂和：

S(n,m)=∑i=0n−1imans=n−S(n,m+1)S(n,m)\begin{aligned}S(n,m)&=\sum_{i=0}^{n-1}i^m\\ans&=n-\dfrac{S(n,m+1)}{S(n,m)}\end{aligned}

考虑如何计算它：观察到 m≤5000m\le5000，可以考虑拉格朗日插值（一个 mm 次幂和的前缀和一定是一个不超过 m+1m+1 次幂的多项式，计算出来前 m+2m+2 即可）：

记：

xi=iyi=S(i,m)k=n\begin{aligned}x_i&=i\\y_i&=S(i,m)\\k&=n\end{aligned}

转化为基本的拉格朗日插值：

已知：ff 是一个次数不超过 nn 的多项式，f(xi)=yi(i∈[0,n]⋂N)f(x_i)=y_i(i\in[0,n]\bigcap\mathbf{N})，求：f(k)=?f(k)=?。

解

f(k)=∑i=0nyi∏i≠jk−xjxi−xjf(k)=\sum^n_{i=0}y_i\prod_{i\ne j}\frac{k-x_j}{x_i-x_j}

其中，展开得：

f(k)=y0(k−x1x0−x1 ⋅ k−x2x0−x2 ⋯k−xnx0−xn )+y1(k−x0x1−x0 ⋅ k−x2x1−x2 ⋯k−xnx1−xn )+y2(k−x0x2−x0 ⋅ k−x1x2−x1 ⋯k−xnx2−xn )⋯+yn(k−x0xn−x0 ⋅ k−x1xn−x1 ⋯k−xn−1xn−xn−1 )\begin{aligned}f(k)&=y_0\left(\frac{k-x_1}{x_0-x_1}~\cdot~\frac{k-x_2}{x_0-x_2}~\cdots\frac{k-x_n}{x_0-x_n}~\right)\\&+y_1\left(\frac{k-x_0}{x_1-x_0}~\cdot~\frac{k-x_2}{x_1-x_2}~\cdots\frac{k-x_n}{x_1-x_n}~\right)\\&+y_2\left(\frac{k-x_0}{x_2-x_0}~\cdot~\frac{k-x_1}{x_2-x_1}~\cdots\frac{k-x_n}{x_2-x_n}~\right)\\&\cdots\\&+y_n\left(\frac{k-x_0}{x_n-x_0}~\cdot~\frac{k-x_1}{x_n-x_1}~\cdots\frac{k-x_{n-1}}{x_n-x_{n-1}}~\right)\end{aligned}

该式对所有 i∈[0,n]⋂Ni\in[0,n]\bigcap\mathbf{N} 均满足 f(xi)=yif(x_i)=y_i。

理由是 yj(j≠i)y_j(j\ne i) 的系数必含 k−xik-x_i 项，此时求 f(k)f(k)：

f(k)=y0(k−xix0−xi ⋯ )+y1(k−xix1−xi ⋯ )+y2(k−xix2−xi ⋯ )⋯+yn(k−xixn−xi ⋯ )\begin{aligned}f(k)&=y_0\left(\frac{k-x_i}{x_0-x_i}~\cdots~\right)\\&+y_1\left(\frac{k-x_i}{x_1-x_i}~\cdots~\right)\\&+y_2\left(\frac{k-x_i}{x_2-x_i}~\cdots~\right)\\&\cdots\\&+y_n\left(\frac{k-x_i}{x_n-x_i}~\cdots\right)\end{aligned}

这些项用 k=xik=x_i 代入即全部消掉，其中只有一项 yiy_i 的系数是特殊的：

f(k)=yi(xi−x0xi−x0 ⋅ xi−x1xi−x1 ⋯ xi−xnxi−xn)f(k)=y_i\left(\frac{x_i-x_0}{x_i-x_0}~\cdot~\frac{x_i-x_1}{x_i-x_1}~\cdots~\frac{x_i-x_n}{x_i-x_n}\right)

后面一项的值显然为 11，因此对任意 i∈[0,n]⋂Ni\in[0,n]\bigcap\mathbf{N} 均有 f(xi)=yif(x_i)=y_i。

又根据 n+1n+1 个点能唯一确定一个不超过 nn 次的多项式，该 f(k)f(k) 的表达式可以表示所有的 f(k)f(k)。

直接将题中 kk 的值代入解析式 f(k)=∑i=0nyi∏i≠jk−xjxi−xjf(k)=\sum^n_{i=0}y_i\prod_{i\ne j}\frac{k-x_j}{x_i-x_j} 求解即可。

考虑证明，严谨证明直接套贝叶斯公式就好，其实没那么麻烦：

在第 ii 名的分布列的比例关系是：0m,1m,⋯ ,(n−1)m0^m, 1^m, \cdots, (n-1)^m，然后就直接对着算就好了。

#include<cstdio>
#include<cstring>
#define int long long
int init(){
	char c = getchar();
	int x = 0, f = 1;
	for (; c < '0' || c > '9'; c = getchar())
		if (c == '-') f = -1;
	for (; c >= '0' && c <= '9'; c = getchar())
		x = (x << 1) + (x << 3) + (c ^ 48);
	return x * f;
}
void print(int x){
	if (x < 0) x = -x, putchar('-');
	if (x > 9) print(x / 10);
	putchar(x % 10 + '0');
}
const int N = (int) 1e4 + 5, Mod = 998244353;
int x[N], y[N]; int ok[N << 1];
int quick_mod(int a, int b){ // 这道题 5000 * 5000 的数据范围卡朴素快速幂，需要优化普通的快速幂到 O (M * log + m^2)，第一个 M 是可能的逆元个数，本题中约为 2m 个
    if (b == Mod - 2) {
        if (a < N && ok[a + N] != -1)
            return ok[a + N];
    }
	int s = 1, olda = a, oldb = b;
    a = (a + Mod) % Mod;
	while (b) {
		if (b & 1) s = s * a % Mod;
		a = a * a % Mod; b >>= 1;
	}
    if (oldb == Mod - 2 && olda < N)
        ok[olda + N] = s;
	return s;
}
int S(int n, int m){
	for (int i = 1; i <= m + 2; ++i)
		x[i] = i, y[i] = (y[i - 1] + quick_mod(i - 1, m)) % Mod;
	int k = n, ans = 0;
	for (int i = 1; i <= m + 2; ++i) {
		int mul = y[i];
		for (int j = 1; j <= m + 2; ++j)
			if (i != j)
                mul = (mul * ((k - x[j]) % Mod + Mod) % Mod) * quick_mod(x[i] - x[j], Mod - 2) % Mod;
		ans = (ans + mul) % Mod;
	}
	return ans;
}
signed main(){
    memset(ok, -1, sizeof(ok));
	int n = init(), m = init();
    int s1 = S(n, m + 1), s2 = S(n, m);
	int ans = s1 * quick_mod(s2, Mod - 2) % Mod;
	print(((n - ans) % Mod + Mod) % Mod), putchar('\n');
}