?fitdistr it says you can pass it a density function (CDF), so we just need to define a density function for the Burr distribution. Note that I'm not familiar with the Burr distribution so I just pulled it's CDF off of Wikipedia.
dburr <- function(x, c = 1, k = 1) 1 - (1 + x ^ c) ^ (-k)
# Simulate data from log logistic for a test case
obs <- rllog(100)
fitdistr(x = obs,
densfun = dburr,
start = list(c = 1, k = 1), # need to provide named list of starting values
lower = list(c = 0, k = 0)) # and named list of lower bounds since c, k > 0
I'm also not sure what the relationship is between Burr and log-logistic, so I don't know what the "right" answer is...
Looking under the CRAN Task View on Distributions, apparently the
VGAM package includes the Pareto Type-IV distribution, which includes Burr's distribution somehow. So if you know how to parameterize Pareto-IV to become Burr, you can use their
dparetoIV function to for
fitdstr, and their
rparetoIV if you want to simulate data.