I cannot wrap my head around the concept of accuracy as a non-differentiable reward function. Do we need to find the function and then check if it is mathematically non-differentiable?
In reinforcement learning (RL), the reward function is often not differentiable with respect to any learnable parameters. In fact it is quite common to not know the function at all, and apply a model-free learning method based purely on sampling many observations of state, action, reward and next state. For reward, you don't need to know the reward function as a function of your parameters, but you do need to know how to caclulate it from each observation. That means that the reward is either provided as part of the environment, or it is clear how to calculate it from the initial state, action or next state. In the case of an agent that creates a neural network then reports the accuracy, you can probably view the NN training process and its results including accuracy on a validation set as a "black box" that reports back an arbitrary reward signal.
In Reinforcement Learning, when reward function is not differentiable, a policy gradient algorithm is used to update the weights of a network
Policy gradient methods are one broad family of RL methods which are also often model-free. Alternatively, value-based methods are broadly the other choice (e.g. Q-learning) or can be used in combination with policy gradient methods (e.g. Actor-Critic). All of these can be model-free, and often are.
Unless you want to apply a model-based RL method, then you have no need to find the reward function in an explicit form based on any parameters. Even if you did want to use a model-based approach, the reward function does not need to be differentiable.
I was wondering if I can use another value, for example silhouette score (in a different scenario) as the reward signal?
Yes this would likely be viable, and could be made to work similarly to the paper, provided there is some meaningful connection between the parameters you manipulate and the end result.
RL can be made to work as a generic way of optimising numeric measure influenced due to a vector of parameters. That includes indirect optimisation of neural networks by non-differentiable metrics such as accuracy. There is nothing special about the choice of accuracy here, the metric just needs to be related to the parameters being learned.
There is a catch. Indirect learning using RL can produce very noisy training data, due to e.g. covariance between all the parameters, and it can take many samples to get a clear gradient signal through all the noise and make meaningful updates. It can be inefficient compared to other optimisation methods if they are available - however, the ability to extract a gradient for the parameters even when rules are complex and non-differentiable is a nice feature of RL.