Running a regression on data that is spatially autocorrelated is fine, and unavoidable in most scenarios (e.g. ecological modelling).
It is when you have SAC in your residuals that you have issues. The assumptions of independence are not met and the chance of Type 1 error is increased. Not to mention potential for unstable/biased parameter estimates.
Use a virtual layer:
Layer -> Create Layer -> New Virtual Layer
and enter a query like this:
SUM(CASE WHEN h.Hyear>2000 THEN 1 ELSE 0 END) "H_New",
SUM(CASE WHEN h.Hyear<=2000 THEN 1 ELSE 0 END) "H_Old",
SUM(CASE WHEN h.Fyear>2000 THEN 1 ELSE 0 END) "F_New",
SUM(CASE WHEN h.Fyear<=2000 THEN 1 ELSE 0 END) "...
Use "Select by Location" under Vector > Research Tools. Use that to select all of the points which are in the polygon.
Then open the attribute table of the points, click "Select features using an expression", use something like the following:
"type" = 'individual' AND "year of construction" < 2000
Then click the small drop down arrow next to Select ...
Before you address the second part of the question, which will likely be solved using raster algebra or a similar function, you will need to clarify in statistical terms what you mean by "probability of occurrence".
The "Kernel Density" tool gives you a raster where the cell values represent the "number of seals per unit area". In other words, how many ...