Your LAD function:

LADf=function(par) {(sum(y-par[1]-par[2]*x)^2)}

Looks exactly like least squares for me. Hence what you do is minimizing sum of squares, not abs of deviance. You need something like

LADf <- function(par) { sum(abs(y - par[1] - par[2]*x)) }

Note that this function is not differentiable, so you have to use an optimizer that can handle that (Nelder-Mead or SANN for instance).

LAD estimator is also equivalent to median regression so you can do this with quantile regression package instead.

The definition of "analogue" is not clear.

If you view $R^2$ as a metric to measure the "goodness of the fit" in regression setting. Then, likelihood can be used to evaluate the "goodness of the fit" for LDA (Linear Discriminant Analysis).

EDIT: I think OP was trying to ask LAD (Least absolute deviations regression) but not LDA. Here is my answer to LAD.

I would suggest OP to review the loss function in regression setting. For squared loss regression, the loss function is

$$\sum_i (y_i-\hat y_i)^2$$

where $\hat y_i$ is the predicted value for data $i$, and the sum is over all the data points.

For least absolute deviation loss, the loss function is

$$\sum_i |(y_i-\hat y_i)|$$

Therefore, to evaluate the "goodness of fit", we can exam the loss value. Now, the problem is with the loss value metric, it will depend on number of data points, i.e., the more data points we have, in general, the higher value the will be in aforementioned 2 formulas.

One way of dealing it is divided by number of data points. And if we want to do one step more, we can "normalize" it into 0 to 1 using RSS, which is R square.

But I think dividing the number of data points should be good for your purpose.

Your attempt didn't work because you specified "qr" as the method for geom_smooth rather than "rq", which gives you a least absolute deviation best fit line.

ggplot(data, aes(x=country_mean_rep, y=country_mean_crime)) + 
  geom_point() + 
  geom_smooth(aes(colour="Linear", fill="Linear"), 
              method="lm", 
              formula=y ~ x,
              se = F) +
  geom_smooth(aes(colour="Least Absolute Deviation", fill="Least Absolute Deviation"), 
              method="rq",
              se = F,
              formula = y ~ x) +
  labs(colour="Functional Form", fill="Functional Form") +
  geom_text(aes(label=country_name), nudge_y=0.02) +
  theme_bw()

Output:

The caveat of this solution is that for it to work, you have to set se = F for your least absolute deviation line, which means you can't plot confidence intervals for this model. You can still have your confidence interval for the linear model if you like; just set se = T in the respective geom_smooth call.