Here's an alternative solution:
FUDGE: {-1, 0, +1}
function: ROLL:s reroll up to SKILL:n {
N: [lowest of SKILL and [count -1 in ROLL]]
result: NdFUDGE + {1 .. #ROLL-N}@ROLL
}
loop SKILL over {0..4} {
output [4dFUDGE reroll up to SKILL] named "skill [SKILL]"
}
The function should be mostly self-explanatory; the one part that may require explanation is {1 .. #ROLL-N}@ROLL
, which sums all but the last N
elements of the sequence ROLL
. By default, AnyDice sorts dice rolls in descending numerical order, so the last elements are the lowest.
In graph mode, the outputs of this program look like this:

Note how the differences between skill levels 2, 3 and 4 are fairly minor, since rolling three or four -1s on 4dF is quite unlikely to begin with.
BTW, the program above assumes, as you say at the end of your question, that players are conservative and will only reroll negative rolls. If your players like to take risks, they might decide to reroll zeros as well, in which case the results would look like this instead:

Note how the averages are still the same, but the results for higher skills have a lot more variance. In particular, the probabilities of rolling a perfect four with a positive skill are a lot higher this way.
(The only difference between the programs used to generate the two graphs above are that the second one uses [count {0, -1} in ROLL]
instead of [count -1 in ROLL]
.)
In particular, if your players are trying to roll against a specific minimum target number, it may make sense for them to only roll as many zeros as needed to maximize their chance of meeting the target.
The optimal strategy in this cases depends on whether the players can reroll the dice one by one, and decide after each roll whether they want to continue rerolling, or whether they have to first decide which dice they want to reroll and then roll them all at once.
In the first case (i.e. sequential rerolls) the optimal decision-making process can be simulated with a recursive AnyDice function:
FUDGE: {-1, 0, +1}
function: first N:n of SEQ:s {
FIRST: {}
loop I over {1..N} { FIRST: {FIRST, I@SEQ} }
result: FIRST
}
function: ROLL:s reroll up to SKILL:n target TARGET:n {
if ROLL + 0 >= TARGET { result: 1 } \- success -\
if #ROLL = 0 | SKILL = 0 | #ROLL@ROLL = 1 { result: 0 } \- failure -\
FIRST: [first #ROLL-1 of ROLL]
result: [FIRST reroll up to SKILL-1 target TARGET - 1dFUDGE] \- reroll -\
}
loop TARGET over {-3..4} {
loop SKILL over {0..4} {
output [4dFUDGE reroll up to SKILL target TARGET] named "target [TARGET], skill [SKILL]"
}
}
Here, the main function ROLL reroll up to SKILL target TARGET
returns 1 if the given roll is equal to or greater than the target, and 0 if it's less than the target and no improvement is possible (i.e. there are no more dice left in the pool, no more rerolls are allowed or the lowest die is already a +1). Otherwise it removes the lowest die from the pool (using a helper function, since AnyDice doesn't happen to have a suitable one built into it), decreases the number of remaining rerolls by one, subtracts 1dF from the target value to simulate a single reroll and then calls itself recursively.
The output of this program is a bit awkward to parse from AnyDice's normal bar / line graph view, so I instead exported it and ran it through the Python script from this earlier answer to turn it into a nice two-dimensional grid that I could import into Google Sheets. The results, as a heat map and as a multi-bar graph, look like this:

In the second case (i.e. all rerolls at once) we first need to figure out what the optimal strategy actually is. A moment's thought shows that:
One should always reroll any -1s, since doing so can never decrease the result. Since the expected average result of a reroll is 0, the expected average after rerolling all -1s equals the number of +1s in the initial roll.
Rerolling a zero doesn't change the expected average result, but it does increase the variance, i.e. it makes the actual result more likely to be further away from the average in either direction. Thus, one should only reroll zeros if the expected average result after rerolling all -1s (i.e. the number of +1s in the initial roll) is below the target number.
Applying this logic in AnyDice results in something like this program:
FUDGE: {-1, 0, +1}
function: ROLL:s reroll up to SKILL:n target TARGET:n {
if [count +1 in ROLL] >= TARGET {
N: [lowest of SKILL and [count -1 in ROLL]]
} else {
N: [lowest of SKILL and [count {0, -1} in ROLL]]
}
result: (NdFUDGE + {1 .. #ROLL-N}@ROLL) >= TARGET
}
loop TARGET over {-3..4} {
loop SKILL over {0..4} {
output [4dFUDGE reroll up to SKILL target TARGET] named "target [TARGET], skill [SKILL]"
}
}
Exporting the output of this script and running it through the same Python script and spreadsheet gives the following heat map and bar graph:

As you can see, the results aren't actually that different from the sequential rerolls case. The biggest differences occur with high skills and intermediate target numbers: for example, with a skill of 4, being able to perform the rerolls one at a time and stop at any point raises the average success rate from 75.3% to 81% for a target of +1, or from 51.6% to 58.3% for a target of +2.
Ps. I did manage to figure out a way to make AnyDice collect the "success rate vs. target" values from the two programs above into a single distribution for each skill value, allowing them to be drawn directly by AnyDice as bar charts or line graphs (in "at least" mode) without having to use Python or spreadsheets.
Unfortunately, the AnyDice code to do that is anything but simple. The hardest(!) part turned out to be finding a way to make AnyDice subtract two probabilities (e.g. 1/2 − 1/3 = 1/6). The best way I know of to perform this seemingly trivial task in AnyDice involves nontrivial manipulation of conditional probabilities and an iterated loop. And it crashes AnyDice if you try to calculate 0 − 0 with it.*
Anyway, just for completeness, here's the AnyDice code for calculating and plotting the distribution of the "highest beatable target" for various skill levels (and for each of the two rerolling mechanics described above) with some comments added for readability:
\- predefine a fudge die -\
FUDGE: d{-1, 0, +1}
\- miscellaneous helper functions used in the code below -\
function: first N:n of SEQ:s {
FIRST: {}
loop I over {1..N} { FIRST: {FIRST, I@SEQ} }
result: FIRST
}
function: exclude RANGE:s from ROLL:n {
if ROLL = RANGE { result: d{} } else { result: ROLL }
}
function: sign of NUM:n {
result: (NUM > 0) - (NUM < 0)
}
function: if COND:n then A:d else B:d {
if COND { result: A } else { result: B }
}
\- a helper function to subtract two probabilities (given as {0,1}-valued dice) -\
function: P:d minus Q:d {
DIFF: P - Q
loop I over {1..20} {
TEMP: [exclude 0 from DIFF]
DIFF: (DIFF != 0) * [sign of TEMP + TEMP]
}
result: [exclude -1 from DIFF]
}
\- this function calculates the probability of meeting or exceeding the target -\
\- value, assuming that each die in the initial roll can be rerolled once and -\
\- that the player may stop rerolling at any point -\
function: ROLL:s reroll one at a time up to SKILL:n target TARGET:n {
if ROLL + 0 >= TARGET { result: 1 } \- success -\
if #ROLL = 0 | SKILL = 0 | #ROLL@ROLL = 1 { result: 0 } \- failure -\
FIRST: [first #ROLL-1 of ROLL] \- remove last (=lowest) original roll -\
TNEW: TARGET - 1dFUDGE \- adjust target value depending on reroll -\
result: [FIRST reroll one at a time up to SKILL-1 target TNEW] \- reroll -\
}
\- this function calculates the probability of meeting or exceeding the target -\
\- value, assuming that each die in the initial roll can be rerolled once but -\
\- the player must decide in advance how many of the dice they'll reroll; the -\
\- optimal(?) decision rule in this case is to always reroll all -1s and to -\
\- also reroll 0s if and only if the number of +1s in the initial roll is less -\
\- than the target number -\
function: ROLL:s reroll all at once up to SKILL:n target TARGET:n {
if [count +1 in ROLL] >= TARGET {
N: [lowest of SKILL and [count -1 in ROLL]]
} else {
N: [lowest of SKILL and [count {0, -1} in ROLL]]
}
result: (NdFUDGE + {1 .. #ROLL-N}@ROLL) >= TARGET
}
\- this function collects the success probabilities given by the two functions -\
\- above into a single custom die D, such that the probability that D >= N is -\
\- equal to the probability of the player meeting or exceeding the target N; -\
\- the SEQUENTIAL flag controls which of the functions above is used -\
function: collect results for SKILL:n from MIN:n to MAX:n sequential SEQUENTIAL:n {
BOGUS: MAX + 1
DIST: 0
PREV: 1
loop TARGET over {MIN..MAX} {
if SEQUENTIAL {
PROB: [4dFUDGE reroll one at a time up to SKILL target TARGET + 1]
} else {
PROB: [4dFUDGE reroll all at once up to SKILL target TARGET + 1]
}
DIST: [if d{MIN..TARGET} < TARGET then DIST else [if [PREV minus PROB] then TARGET else BOGUS]]
PREV: PROB
}
result: [exclude BOGUS from DIST]
}
\- finally we just loop over possible skill values and output the results -\
loop SKILL over {0..4} {
output [collect results for SKILL from -4 to 4 sequential 1] named "skill [SKILL], one at a time"
}
loop SKILL over {0..4} {
output [collect results for SKILL from -4 to 4 sequential 0] named "skill [SKILL], all at once"
}
and a screenshot of the output (in "at least" line graph mode):

A note on the interpretation of the output generated by the program above: The probability distributions shown on the graph above don't correspond to the results of any single dice rolling strategy; rather, they are artificially constructed distributions (i.e. "custom dice" in AnyDice jargon) such that the probability of rolling at least \$N\$ on a single roll of the custom die equals the probability of the player being able to roll at least \$N\$ on 4dF with the given rerolling mechanic (one at a time vs. all at once) and the given maximum number of rerolls, assuming that the player uses the optimal rerolling strategy for that particular target \$N\$.
In other words, looking at the output in "at least" mode, we can see that a player with skill level 4 has a 51.62% chance of successfully rolling +2 or more (using the all-at-once rerolling mechanic) if they're using their available rerolls in the way that maximizes that particular chance. The output also correctly shows that the same player has a 75.28% chance of rolling +1 or more if they choose to optimize for that instead, but they'll need different rerolling strategies to reach those two goals.
And the "probability" of 23.65% for rolling exactly +1 on the custom die described above really has no sensible meaning, except that it's (approximately, due to rounding) the difference between 75.28% and 51.62%. Which I guess is why it's so hard to calculate with AnyDice. :P I suppose you could interpret it as a measure of how much harder a target of +2 is to meet using the given skill and rerolling mechanic than a target of +1, in some sense, but that's about it.
*) That crash might be related to what I'm pretty sure is a bug in AnyDice that I found while developing this code, causing one of my early test programs to generate really weird output with things like 97284.21% probabilities(!). The test program also eventually crashes if you increase the iteration count further.