Assume I have two datasets, each one containing 5000 samples, and each sample has three dimensions. I am looking for a way to "reorder" the samples in one (or probably both) dataset such that the sum of least square error between these two data sets be minimum.
Please see the attached photo for clarification.
What you have is called as assignment problem (or a marriage problem, you want to marry off your best matches!). It is an instance of combinatorial optimization, and can be solved by the Hungarian algorithm, but in this wiki links you can also find other algorithms.
Below an illustration using R. I sample uniformly 9 points on a circle with radius 1 (red), and other 9 points on a circle with radius 2 (blue), and then we marry off red and blue:
The R code used:
library(RcppHungarian)
N <- 9 ### Problem size
### First a very simple example: N random points on a circle with radius 1,
### then radius two
set.seed(2) # Seed chosen for nice plot!°
angle1 <- runif(N, 0, 2*pi)
angle2 <- runif(N, 0, 2*pi)
RED <- cbind( cos(angle1), sin(angle1) )
BLUE <- cbind( 2*cos(angle2), 2*sin(angle2) )
### Calculating cost (Euclidean distance) matrix:
Cost <- matrix(as.numeric(NA), N, N)
for (i in 1:N)
for (j in 1:N) {
Cost[i, j] <- sqrt(sum( (RED[i]-BLUE[j])^2))
}
Sol <- HungarianSolver( Cost )
### The pairings in Sol$pairs
### Plotting:
BLUE_reordered <- BLUE[Sol$pairs[, 2], ]
ARCS <- matrix(as.numeric(NA), 3*N-1, 2)
ARCS[seq(from=1, by=3, length.out=N), ] <- RED
ARCS[seq(from=2, by=3, length.out=N), ] <- BLUE_reordered
plot(BLUE, col="blue", cex=2, pch=19, xlab="", ylab="")
points(RED, col="red", cex=2, pch=19)
lines(ARCS, col="orange")
(I did not try this code with N=5000, as you have, you will need a computer with much memory).