**The setting**
I have a Geopackage point layer in QGIS 3.18 with over 350.000 features for an area of about 1700 km² (extent ca. 50\*60 km), representing centroids of buildings. The points contain an attribute with the construction year of the building: from 1000 to 2020. A few statistical values, based on `Basic statistics for fields`, can be found below. CRS is `EPSG:2056` (projected CRS for Switzerland, units=m).

**What I want to do**
The idea now is to find for each building the *nearest building that is older* and create an attribute `nearest_older` with the `fid` of this next older building - see the visualization at the bottom.

In a conceptual sense, it is similar to the concept of [Topographic isolation](https://en.wikipedia.org/wiki/Topographic_isolation): for a summit, find the minimum distance to a point of equal/higher elevation.

**What I tried**
I found an expression that works in principle, see it below. For some additional ideas how to proceed, see at the bottom as well.

**The problem**
The expression will not be efficient for such a huge data set (at least not in reasonable time). For testing, I reduced the computation time by: A) let it run only for the first 999 features (`if (fid< 1000`) and B) limited the number of nearest features to take into consideration to 1000 (`limit:=1000`). Even with this limitations, calculation runs for about 15 minutes.

**The Question**
How could the expression be improved for performance (calculation time)?
There will be probably no ideal solution, just an optimal one. So I'm interested in a solution that significantly reduces the calclulation time in a way that it becomes operationable for my dataset.

**The expression**

```
if (
fid< 1000,
with_variable (
    'cy',
    "construction_year" ,
    array_first (
        array_filter (  
            overlay_nearest( 
                @layer, 
                $id,  
                limit:=1000
            ),
            attribute (
                get_feature_by_id (
                    @layer, 
                    @element
                ), 
                'construction_year'
            ) < @cy
        )
    )
),
'')

```

**Explanation how the expression works**

1. `overlay_nearest( )` get an ordered array of the `id` for the x nearest neighbours (x here is set to 1000)
2. `array_filter ()`: filter this array in a way that only elements (`@element`) are retained where the `attribute ()` named `construction_year` is smaller than the variable `@cy` - and this variable is defined as the construction year of the current feature. So the filter returns for each feature an array of id's, where the construction year is older than the one of the current feature.
3. `array_first ()` returns the first id of this value. As `overlay_nearest` from step 1 orders the array based on closeness to the current feature, this is the id of the `nearest_older` I am looking for.

**Further ideas how to improve performance**

Here are a few ideas I had how to make the expression more efficient. All of these approaches have some shortcomings, however. And particularly, I'm stuck how to combine them in the most efficient way. And probably there are other approaches, not considered yet.

1. Run it sequentially, like in the if-clause above: for features 1 to 999, than for 999 to 1999 and so on. Not very efficient, however.
2. `limit:=1000` to reduce the number of elements in the array created by `overlay_nearest()` to 1000 is unflexible: the newer a building is, the higher chances are that very close to it, you'll find one that is older. Thus, for the majority of the buildings (who are constructed in the last few decades), you won't need to identify a fixed number of 1000 nearest neighbors - a number of 50 or 100 or so would be OK. So the fixed value `1000` could be replaced by a formula that returns an inverse proportional value regarding the construction year. However, how to get an "optimum" formula, based on the distribution of the values in my field `construction_year`? For this, compare the statistical values below.
3. Not for every feature a "match" has to be found in one pass. For some features, in the first round and using the limits defined, no matching `id` with an older building could be found. These NULL value cases could be calculated (based on a condition `if NULL`) in a next iteration. So an iterative approach could be used - but how to set it up?

**A few statistical values**

These values are based on `Basic statistics for fields` for the field `construction_year`

```
Count: 336856
Unique values: 629
Minimum: 1000
Maximum: 2020
Range: 1020
Mean: 1954
Median: 1971
Standard deviation: 64.5
Coefficient of Variation:  0.033
Minority (rarest occurring value): 1004
Majority (most frequently occurring value): 1980
First quartile: 1935
Third quartile: 1994
Interquartile Range (IQR): 59

```

*Histogram, created with DataPlotly plugin: distribution of construction years, with the highest peak at 1980/81, lower peaks at 1928-1931 etc.:*
[![enter image description here](https://i.stack.imgur.com/BA5Hj.png)](https://i.stack.imgur.com/BA5Hj.png)

**Visualization**

This shows the principle: each point is labeled with its construction year and connected by a red arrow to the nearest point with an older construction year. The point labeled with `1958` at the very bottom is connected to a point with label `1940`, even though it has four neighboring points at a nearer distance, but with newer construction date: `1986`, `1969`, `1996` and `1960` - so it goes on until the first (nearest) point is found with an older construction date:
[![enter image description here](https://i.stack.imgur.com/nJxHy.png)](https://i.stack.imgur.com/nJxHy.png)