I think you may be conflating likelihood with maximum likelihood methods. I'll try to separate them as best I can below.
The likelihood function is one which relates the probability of an observation with a particular parameterization of a specific distributional family. It is not necessarily the same thing as the probability of the observation, as in the continuous case, the likelihood is still well defined at point but the probability is the integral of the pdf over a span and is 0 at every single point. Nevertheless, given both a distributional form (e.g. exponential vs. gamma) and a specific parameterization (a $\beta$ for the exponential or an $\alpha, \beta$ for the gamma) it is a measure of relative probability of one observation versus another.
Maximum Likelihood Estimation
Knowing that, if one has a set of observations and chooses a distributional form, then a logical approach to finding estimates of the parameters of that specific distributional form is to find parameters for which the likelihood of the set of observations is greatest. This is the principle of maximum likelihood estimation (MLE). The actual optimization is almost always performed on the negative log-likelihood for performance and over/underflow reasons, but the principle remains the same. The parameter(s) "most likely" to have generated this data are the ones for which the combined likelihood of all the observations (product of actual or sum of the logs) is the largest.
In MLE, one needs to select a distributional form first and then solve for the parameters. When dealing with count data in particular, often the data exhibits properties that don't allow it to be cleanly modeled by a known distribution with a closed form. For example, the Poisson distribution requires that the variance equal the mean. The negative binomial distribution requires (in most parameterizations at least) that the variance be a function of the square of the mean. What happens when the variance is a linear function of the mean? There is no simple discrete distribution with this property, so how can we find a likelihood without a distributional form? Roger Wedderburn, one of the developers of the GLM framework and the use of IWLS for MLE in that context, also proved that one can use a function similar to a true likelihood—a quasi-likelihood function—despite it not being a "true" distribution so long as the mean-variance relationship is well defined (among other requirements). This allows for using the existing mechanics of IWLS/MLE to be used to create models in the presence of overdispersion and which would otherwise not fit cleanly into the forms of known distributions.
Restricted Maximum Likelihood (REML)
MLE does not always return an unbiased estimator. For example, with the classic normal distribution, the MLE of the mean is the sample mean, which is unbiased, but the MLE of the variance is the population variance, which we know is biased. When using the GLM framework (for which the IWLS approach is actually MLE in disguise) to solve for the variance, often an unbiased estimator is desired. Applying MLE not to the raw data but to a transformation of the raw data, can result in an unbiased estimate. This underlies the use of restricted maximum likelihood. As the transformed data does not necessarily encompass all the data, the method is called "restricted".
The "true" likelihood of a distribution may involve very complicated normalizing factors, especially in multivariate cases. This may make using true maximum likelihood estimation intractable. However, a simplified function of the observations—or a subset of the observations—may be mathematically tractable and allow for estimating a "good" optimum even though it may not be the "best". So while these functions are not the true likelihood, we may treat them as such for the purpose of fitting the distributions. For more information, please see Besag (1975) or Arnold & Strauss (1991).