This problem can be solved in polynomial time. Therefore, it's likely not NP-hard (unless P = NP).
After a little bit of algebraic manipulation, we find that the condition on $S$ is equivalent to requiring that
$$2 P_i \ge \varepsilon + \sum_{j \in S} P_j \quad \text{for all $i \in S$}.$$
This, in turn, is equivalent to requiring
$$2 \min \{P_i : i \in S\} \ge \varepsilon + \sum_{j \in S} P_j.$$
For convenience, sort the $P_i$'s into increasing order, so $P_1 \le P_2 \le \dots \le P_n$.
Now one can prove that if there is any set of cardinality $k$ that satisfies this equation, then there's a set of the form $S = \{t,t+1,t+2,\dots,t+k-1\}$ that satisfies this equation, for some $t$. (This can be proven using an exchange argument, swapping pairs of indices as needed and showing that the swap leaves the left-hand side unchanged and only decreases the right-hand side.)
Now you can try all possible values of $t$ and see whether the corresponding set satisfies the condition. There are only $n$ possible values of $t$, so this can certainly be tested in polynomial time. In fact, with clever programming (prefix sums, etc.), you can implement this in $O(n \lg n)$ time, or in $O(n)$ time if the $P_i$'s are already provided in sorted order.
More generally, if you want to find the set $S$ of maximum cardinality that satisfies the equation, such a set will have the form $S=\{t,t+1,t+2,t+3,\dots,u\}$ for some $t,u$. Therefore, you can try all such possibilities and see which has largest cardinality. This too can be done in $O(n \lg n)$ time using appropriate methods (binary search, prefix sums). Alternatively, you can do binary search on $k$, and apply the decision algorithm of the previous paragraph for each candidate value of $k$ to find the largest $k$ for which a solution exists. The running time is again $O(n \lg n)$.