{"action":"created","issue":{"url":"https://api.github.com/repos/scipy/scipy/issues/4523","labels_url":"https://api.github.com/repos/scipy/scipy/issues/4523/labels{/name}","comments_url":"https://api.github.com/repos/scipy/scipy/issues/4523/comments","events_url":"https://api.github.com/repos/scipy/scipy/issues/4523/events","html_url":"https://github.com/scipy/scipy/issues/4523","id":57528255,"number":4523,"title":"rbf with gaussian kernel seems to produce more noise than original data","user":{"login":"eyurtsev","id":3205522,"avatar_url":"https://avatars.githubusercontent.com/u/3205522?v=3","gravatar_id":"","url":"https://api.github.com/users/eyurtsev","html_url":"https://github.com/eyurtsev","followers_url":"https://api.github.com/users/eyurtsev/followers","following_url":"https://api.github.com/users/eyurtsev/following{/other_user}","gists_url":"https://api.github.com/users/eyurtsev/gists{/gist_id}","starred_url":"https://api.github.com/users/eyurtsev/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/eyurtsev/subscriptions","organizations_url":"https://api.github.com/users/eyurtsev/orgs","repos_url":"https://api.github.com/users/eyurtsev/repos","events_url":"https://api.github.com/users/eyurtsev/events{/privacy}","received_events_url":"https://api.github.com/users/eyurtsev/received_events","type":"User","site_admin":false},"labels":[],"state":"open","locked":false,"assignee":null,"milestone":null,"comments":10,"created_at":"2015-02-12T22:16:42Z","updated_at":"2015-02-16T00:05:10Z","closed_at":null,"body":"This is probably user error, but I can't get rbf interpolation working in scipy. (I tried reading through the code, but it's a bit difficult to figure out what's going on...)\r\n\r\nBasically, the more x data points I add the noisier the interpolation gets:\r\n\r\n![illustration](https://cloud.githubusercontent.com/assets/3205522/6178155/63c36dd4-b2da-11e4-8feb-2f770283abcb.png)\r\n\r\nBelow is a short script showing the problem.\r\n\r\nI used sklearn's KneighborsRegressor to show what I expected to see as my output. I also re-implemented the same thing myself which gives identical output to that of sklearn (not included below.)\r\n\r\nI am on scipy 0.15.1.\r\n\r\n```python\r\n\r\n%pylab inline\r\nimport scipy\r\nfrom scipy.interpolate import Rbf\r\nfrom sklearn.neighbors import KNeighborsRegressor\r\n\r\n# Generate some data\r\nx = np.linspace(0, 10, 50)\r\nz = x + 4.0 * randn(len(x))\r\n\r\n# plot data\r\nscatter(x, z, color='black', label='data')\r\n\r\n# RBF regression using sklearn. Use all neighbors and gaussian kernel\r\ndef weights_func(dist_array):\r\n    return exp(-dist_array**2)\r\n\r\nest = KNeighborsRegressor(n_neighbors=50, weights=weights_func)\r\nest.fit(x.reshape(-1, 1), z)\r\nym = est.predict(x.reshape(-1, 1))\r\nplot(x, ym, color='purple', linestyle='-', linewidth=2, label='rbf using KneighborsRegressor')\r\n\r\n\r\n# Rbf from scipy\r\nrbf = Rbf(x, z, function='gaussian', epsilon=1)\r\nplot(x, rbf(x), label='scipy-rbf')\r\nlegend(loc='center left', bbox_to_anchor=(1.05, 0.5), frameon=False)\r\n```"},"comment":{"url":"https://api.github.com/repos/scipy/scipy/issues/comments/74445097","html_url":"https://github.com/scipy/scipy/issues/4523#issuecomment-74445097","issue_url":"https://api.github.com/repos/scipy/scipy/issues/4523","id":74445097,"user":{"login":"argriffing","id":101416,"avatar_url":"https://avatars.githubusercontent.com/u/101416?v=3","gravatar_id":"","url":"https://api.github.com/users/argriffing","html_url":"https://github.com/argriffing","followers_url":"https://api.github.com/users/argriffing/followers","following_url":"https://api.github.com/users/argriffing/following{/other_user}","gists_url":"https://api.github.com/users/argriffing/gists{/gist_id}","starred_url":"https://api.github.com/users/argriffing/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/argriffing/subscriptions","organizations_url":"https://api.github.com/users/argriffing/orgs","repos_url":"https://api.github.com/users/argriffing/repos","events_url":"https://api.github.com/users/argriffing/events{/privacy}","received_events_url":"https://api.github.com/users/argriffing/received_events","type":"User","site_admin":false},"created_at":"2015-02-16T00:05:10Z","updated_at":"2015-02-16T00:05:10Z","body":"FWIW here's a simpler implementation for the 1d case with epsilon=1, no smoothing.  There are stability issues as mentioned in http://duhctaep.blogspot.com/2010/05/kriging-with-numpy.html.\r\n\r\n```python\r\nfrom __future__ import print_function, division\r\n\r\nfrom functools import partial\r\n\r\nimport numpy as np\r\nimport scipy.interpolate\r\nimport matplotlib.pyplot as plt\r\n\r\ndef distances_1d(a, b):\r\n    return np.absolute(a[np.newaxis, :] - b[:, np.newaxis])\r\n\r\ndef _solve_krige(x, w, u):\r\n    r = distances_1d(x, u)\r\n    B = np.exp(-r*r)\r\n    v = np.dot(B, w)\r\n    return v\r\n\r\ndef krige(x, z):\r\n    r = distances_1d(x, x)\r\n    A = np.exp(-r*r)\r\n    w = scipy.linalg.solve(A, z)\r\n    return partial(_solve_krige, x, w)\r\n\r\ndef main():\r\n    # The solver breaks for this input.\r\n    x = np.linspace(0, 10, 15)\r\n    z = x + 4.0 * np.random.randn(len(x))\r\n\r\n    # plot data\r\n    plt.scatter(x, z, color='black', label='data')\r\n\r\n    u = np.linspace(0, 10, 50)\r\n\r\n    # Rbf from scipy\r\n    rbf = scipy.interpolate.Rbf(x, z, function='gaussian', epsilon=1)\r\n    plt.plot(u, rbf(u), label='scipy-rbf')\r\n\r\n    # custom function that is the same as scipy Rbf for 1d\r\n    f = krige(x, z)\r\n    plt.plot(u, f(u), color='purple', linestyle='-', linewidth=2, label='krige')\r\n\r\n    plt.legend(loc='center left', bbox_to_anchor=(1.05, 0.5), frameon=False)\r\n    plt.show()\r\n\r\nmain()\r\n```"}}