Why does R2 increase with fewer samples using adonis? I've noticed in a large ecological dataset that when I subsample my data, the R2 'proportion of variance' explained by my categorical grouping (output from adonis() within the vegan package) gets larger. That is, when I have fewer samples, that R2 increases. Naively, I think that because more samples equates to having greater statistical power, the idea that more samples underestimates (or fewer samples overestimates) the proportion of variance explained by my variable of interest seems paradoxical to me. Could someone help explain this? Also, does this mean comparing R2 for datasets of different sizes is inappropriate? F statistics go down with fewer numbers, so at least that seems to make sense.
 
Ex.:
birds <- read.csv('https://raw.githubusercontent.com/collnell/lab-demo/master/bird_by_fg.csv')
bird.matrix <- as.matrix(birds[,3:9])
bird.dist <- vegdist(bird.matrix, method='bray')
adonis(bird.dist~DIVERSITY, data=birds, permutations = 999, method="bray")

The full dataset is 32 rows and R2 for that is calculated above. However, the mean R2 scores get larger when subsampling down to 20 (and further down to 15):
set.seed(7*11*13)
R2scores <- vector(mode="numeric")
for (i in 1:50)
{
sampsize=20
randsamp <- sample(1:32, sampsize, replace=F)
bird.matrix <- as.matrix(birds[randsamp,3:9])
bird.dist <- vegdist(bird.matrix, method='bray')
ad1 <- adonis(bird.dist~DIVERSITY, data=birds[randsamp,], permutations = 10000, method="bray")
R2scores <- c(R2scores,ad1$aov.tab[1,5])
}
mean(R2scores)

This is also very evident in my dataset of ~1400 samples. Any help would be greatly appreciated!
 A: You didn't really give enough information in your post for a full answer, and your second code block have an error that may have distorted your results (The index i was used twice. When an index is used to count the loops in a for statement, it shouldn't be assigned within the loop! Also, you should have set the random seed to make the code reproducible, and presented output.) Below is that code corrected:
set.seed(7*11*13) 
R2scores <- vector(mode="numeric")
for(i in 1:100){
size <- 20
randsamp <- sample(1:32, size, replace=FALSE)
bird.matrix <- as.matrix(birds[randsamp,3:9])
bird.dist <- vegdist(bird.matrix, method='bray')
ad1 <- adonis(bird.dist ~ DIVERSITY, data=birds[randsamp,], permutations = 10000, method="bray")
R2scores <- c(R2scores, ad1$aov.tab[1,5])
}

mean(R2scores)

In addition to print the mean, you should also look at the histogram. You will then see a rather large spread, but you are right that the mean value is increasing. I have no experience with adonis, but the probable explanation is overfitting.  With a smaller subsample, it is easier to find a really good fit (that will not generalize well to a new sample.)   
A: I've never used adonis, but in most cases R2 is higher with fewer samples, because the model can overfit - see e.g. https://statisticsbyjim.com/regression/r-squared-too-high/ for more details.
