sklearn gives different trees for scaled and unscaled data

I am not very convinced.  I used the same model with exact same hyperparameters and setting the same random state for a synthetically generated dataset for three scenarios. Firstly, made the tree with unscaled data, then using maxmin scaling and then scaling each input arbitrarily to close to double digits. The tree output was the same for the maxmin scaled and arbitrarily scaled data but different when it comes to the unscaled data. I am enclosing the code used for reference. This issue is not observed when I tried the same approach with other machine learning packages like MATLAB. I am also enclosing the zip file with the jupyter notebook and the output trees for reference.
[test_tree.zip](https://github.com/scikit-learn/scikit-learn/files/10844429/test_tree.zip)

```python
import pandas as pd
import numpy as np
from sklearn import tree
import dtreeviz

dfe  = np.random.exponential(scale=1e-7, size=(100000,1))
mld  = np.random.exponential(scale=1e2, size=(100000,1))
nh4  = np.random.exponential(scale=1e-5, size=(100000,1))
no3  = np.random.exponential(scale=1e-5, size=(100000,1))
po4  = np.random.exponential(scale=1e-5, size=(100000,1))
rsn  = np.random.exponential(scale=1e2, size=(100000,1))
si   = np.random.exponential(scale=1e-4, size=(100000,1))
sos  = np.random.exponential(scale=1e2, size=(100000,1))
tos  = np.random.exponential(scale=1e2, size=(100000,1))
wo5  = np.random.exponential(scale=1e-7, size=(100000,1))
phyc = np.random.exponential(scale=1e-4, size=(100000,1))

columns=['dfe', 'mld', 'nh4', 'no3', 'po4', 'rsn','si',
            'sos', 'tos', 'wo5','phyc']
arr = np.concatenate((dfe,mld,nh4,no3,po4,rsn,si,sos,tos,wo5,phyc), axis=1)
da = pd.DataFrame(arr, columns=columns)
df = da.copy(deep=True)

from sklearn.model_selection import train_test_split
from sklearn.metrics import mean_squared_error

X = df.drop('phyc', axis=1)
y = df['phyc']
X_train, X_test, y_train, y_test = train_test_split(X, y,test_size=0.1,random_state=101,shuffle=False)

import graphviz
from sklearn.tree import DecisionTreeRegressor
model = DecisionTreeRegressor(
    criterion='squared_error',
    splitter='best',
    min_samples_leaf=5,
    max_depth=5,
    # min_samples_split=2,
    random_state=101
    )

#Unscaled Tree

unscaled_tree=model.fit(X_train,y_train)
y_pred = unscaled_tree.predict(X_test)
rmse_unscaled = mean_squared_error(y_test,y_pred, squared=False)
print('rmse :',rmse_unscaled)
print(rmse_unscaled /df['phyc'].std())

from sklearn.tree import export_graphviz
estimator=unscaled_tree
export_graphviz(estimator,
                out_file='test_unscaled.dot',
                feature_names=X.columns,
                class_names='phyc',
                rounded=True,proportion=False,
                precision=2, filled=True,)

from subprocess import call
call(['dot','-Tpng', 'test_unscaled.dot', '-o', 'test_unscaled.png','-Gdpi=600'])


from IPython.display import Image
Image(filename = 'test_unscaled.png')

#Tree maxmin scaled

from sklearn.preprocessing import MinMaxScaler
scaler = MinMaxScaler(feature_range=(0,1))
dfs = scaler.fit_transform(da)
dfs = pd.DataFrame(dfs, columns=['dfe', 'mld', 'nh4', 'no3', 'po4', 'rsn','si','sos', 'tos', 'wo5','phyc'])
Xs = dfs.drop(['phyc'], axis=1)
ys = dfs['phyc']

Xs_train, Xs_test, ys_train, ys_test = train_test_split(Xs, ys, shuffle=False,test_size=0.1, random_state=101)

scaled_tree= model.fit(Xs_train,ys_train)
ys_pred = scaled_tree.predict(Xs_test)

from sklearn.tree import export_graphviz
estimator=scaled_tree
export_graphviz(estimator,
                out_file='test_scaled.dot',
                feature_names=X.columns,
                class_names='phyc',
                rounded=True,proportion=False,
                precision=2, filled=True,)
from subprocess import call
call(['dot','-Tpng', 'test_scaled.dot', '-o', 'test_scaled.png','-Gdpi=600'])
from IPython.display import Image
Image(filename = 'test_scaled.png')

#Tree arbitrary scaled

dfr = da.copy(deep=True)
dfr['dfe']=dfr['dfe']*1e9
dfr['nh4']=dfr['nh4']*1e7
dfr['no3']=dfr['no3']*1e7
dfr['phyc']=dfr['phyc']*1e5
dfr['po4']=dfr['po4']*1e6
dfr['si']=dfr['si']*1e5
dfr['wo5']=dfr['wo5']*1e9

Xr = dfr.drop('phyc', axis=1)
yr = dfr['phyc']
Xr_train, Xr_test, yr_train, yr_test = train_test_split(Xr, yr,test_size=0.1,random_state=101,shuffle=False)

arscaled_tree=model.fit(Xr_train,yr_train)
yr_pred = arscaled_tree.predict(Xr_test)

from sklearn.tree import export_graphviz
estimator=arscaled_tree

export_graphviz(estimator,
                out_file='test_arscaled.dot',
                feature_names=X.columns,
                class_names='phyc',
                rounded=True,proportion=False,
                precision=2, filled=True,)


from subprocess import call
call(['dot','-Tpng', 'test_arscaled.dot', '-o', 'test_arscaled.png','-Gdpi=600'])

from IPython.display import Image
Image(filename = 'test_arscaled.png')
```

_Originally posted by @sandupal in https://github.com/scikit-learn/scikit-learn/issues/25495#issuecomment-1447103075_

COMMENT:
We would like to know whether sort algorithm used to construct the tree produces a unique sort or whether it is approximate (i.e. is there some tolerance that is applied). Again, *tree methods should not care about the scale of the inputs* and produce fundamentally different results. We find that when applied to identical datasets scikit-learn does produce such different results but Matlab does not. So it cannot be as simple as "rounding".  We believe it is worth understanding why this is the case- even though we have figured out how to work around it others will not immediate realize that is an issue.

Labels: Needs Triage