Comparing Kendall's Tau values I have a single ordinal outcome called RATING, and a series of continuous covariates (P1, P2, P3, etc). I have used Kendall's Tau statistic to assess the concordance/correlation between RATING and each covariate as per below.
#load Kendall package
library(Kendall)

#create example data
RATING<-sample(1:4,200,replace=TRUE)
P1<-RATING+rnorm(200,0,2.5)
P2<-RATING+rnorm(200,0,3)
P3<-RATING+rnorm(200,0,3.5)

#Calculate Kendall's Tau for each predictor
Kendall(RATING,P1)
Kendall(RATING,P2)
Kendall(RATING,P3)

The resulting Tau statistics will differ in size and significance obviously. My question though is whether there a way to determine which are significantly different from others? For example, can I use the results of these tests to say P1 is significantly more strongly correlated with RATING than P3, or that P3 is significantly less correlated with RATING than all other covariates?
I've found discussion of a method to compare Pearson's R values here How to compare the strength of two Pearson correlations? but I believe it assumes the compared R values are from independent samples, its not clear if it can be applied to Kendall's tau values, and its not clear whether its suitable for multiple comparisons.
Any suggestions greatly appreciated
 A: Here's a rough and ready approach using bootstrapping, borrowing a method from Guillame Rousselet. 
First I'll replicate your data, modifying P1 slightly such that it will have a strong relationship with the Rating just to illustrate. I'm using the the cor function from the default stats package with the method = "kendall" argument throughout.
set.seed(6000)
Rating <- sample(1:4, 200, replace=TRUE)
df <- data.frame( Rating,
                  P1 = Rating + rnorm(200, 0, .8),
                  P2 = Rating + rnorm(200, 0, 3.5),
                  P3 = Rating + rnorm(200, 0, 3))

#Check the 3 Kendall Tau values
sapply(2:4, function(x) cor(df[ , 1], df[ , x], method = "kendall"))
# [1] 0.6535829 0.1738252 0.2347800

So far so good. 
Now we'll want to iteratively resample observations in the dataset, calculate the three overlapping Kendall's $\tau$ values, and take the differences between each pair ($\tau_i-\tau_j)$. Determining whether the bootstrapped distributions of these differences overlap zero provides your test. 
library(boot)
library(broom)

#Create Bootstrap Function
compKendall <- function(data, indices){
  df <- data[indices, ] #Resample over rows

  #Calculate all 3 overlapping correlations
  ks <- sapply(2:4, function(x) cor(df[ , 1], df[ , x], method = "kendall"))

  #Arrange pairwise combinations of correlations
  kPairs <- combn(ks, 2)

  #Return the differences
  return( kPairs[1 , ] - kPairs[2 , ] )
}

#Run the bootstrap with 1000 resamples
kDiff <- boot(df, statistic = compKendall, R = 1000)

#Use broom to get a neat dataframe of differences w/ 95% CI's
broom::tidy(kDiff, conf.int = TRUE, conf.method = "perc")

##  A tibble: 3 x 5
#   statistic     bias std.error conf.low conf.high
#      <dbl>    <dbl>     <dbl>    <dbl>     <dbl>
# 1    0.480   0.00102    0.0533    0.378    0.591 
# 2    0.419  -0.00308    0.0525    0.312    0.518 
# 3   -0.0610 -0.00410    0.0694   -0.205    0.0720

The dataframe shows the estimated differences between Kendall's $\tau$ values, organised as $\tau$ number 1 - 2, 1 - 3, and 2 - 3. 
As expected, the first correlation ($\tau_1 = .65$) appears larger than the other two - CI's for the differences exceed zero. However we couldn't conclude there is a significant difference between the second ($\tau_2 = .17$) and third ($\tau_3 = .23$) overlapping correlations.
