{"action":"created","issue":{"url":"https://api.github.com/repos/sherpa/sherpa/issues/184","repository_url":"https://api.github.com/repos/sherpa/sherpa","labels_url":"https://api.github.com/repos/sherpa/sherpa/issues/184/labels{/name}","comments_url":"https://api.github.com/repos/sherpa/sherpa/issues/184/comments","events_url":"https://api.github.com/repos/sherpa/sherpa/issues/184/events","html_url":"https://github.com/sherpa/sherpa/issues/184","id":136072080,"number":184,"title":"Low-level API example for user stat / dataset","user":{"login":"cdeil","id":852409,"avatar_url":"https://avatars.githubusercontent.com/u/852409?v=3","gravatar_id":"","url":"https://api.github.com/users/cdeil","html_url":"https://github.com/cdeil","followers_url":"https://api.github.com/users/cdeil/followers","following_url":"https://api.github.com/users/cdeil/following{/other_user}","gists_url":"https://api.github.com/users/cdeil/gists{/gist_id}","starred_url":"https://api.github.com/users/cdeil/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/cdeil/subscriptions","organizations_url":"https://api.github.com/users/cdeil/orgs","repos_url":"https://api.github.com/users/cdeil/repos","events_url":"https://api.github.com/users/cdeil/events{/privacy}","received_events_url":"https://api.github.com/users/cdeil/received_events","type":"User","site_admin":false},"labels":[],"state":"open","locked":false,"assignee":null,"milestone":null,"comments":5,"created_at":"2016-02-24T14:10:07Z","updated_at":"2016-03-03T11:34:14Z","closed_at":null,"body":"Is there an example script or IPython notebook that shows how to use the low-level, class-based API to implement a user statistic and/or user dataset as part of a joint fit with other datasets?\r\n\r\nWe want to do a joint likelihood fit for multi-mission datasets (e.g. flux points or OGIP data) that can be represented as Sherpa datasets, as well as some datasets (e.g. Fermi-LAT) for which we'd like to evaluate the likelihood with another software (e.g. the Fermi-LAT Science tools).\r\n\r\nTo make it concrete, let's say we have some flux points where the built-in Data1D and Chi2DataVar are appropriate as in this example: http://python4astronomers.github.io/fitting/low-level.html\r\n\r\nAnd as a second dataset, we have Fermi data and would like to do a joint fit of a Sherpa spectral model to both datasets. For Fermi, the interface should basically be a user-defined stat function we write that is called with a tabulated spectrum (i.e. `energy` and `flux` 1-dim arrays with flux evaluated by Sherpa) and the output is the likelihood on the Fermi dataset.\r\n\r\nIs it clear?\r\n\r\nCan you point us to the relevant classes or even make a toy example to get us started?\r\n\r\nFor context of this question, see discussion here:\r\nhttps://github.com/open-gamma-ray-astro/gamma-astro-data-formats/issues/4#issuecomment-186715011\r\n\r\ncc @giacomov  @tibaldo"},"comment":{"url":"https://api.github.com/repos/sherpa/sherpa/issues/comments/191719128","html_url":"https://github.com/sherpa/sherpa/issues/184#issuecomment-191719128","issue_url":"https://api.github.com/repos/sherpa/sherpa/issues/184","id":191719128,"user":{"login":"dtnguyen2","id":12256848,"avatar_url":"https://avatars.githubusercontent.com/u/12256848?v=3","gravatar_id":"","url":"https://api.github.com/users/dtnguyen2","html_url":"https://github.com/dtnguyen2","followers_url":"https://api.github.com/users/dtnguyen2/followers","following_url":"https://api.github.com/users/dtnguyen2/following{/other_user}","gists_url":"https://api.github.com/users/dtnguyen2/gists{/gist_id}","starred_url":"https://api.github.com/users/dtnguyen2/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/dtnguyen2/subscriptions","organizations_url":"https://api.github.com/users/dtnguyen2/orgs","repos_url":"https://api.github.com/users/dtnguyen2/repos","events_url":"https://api.github.com/users/dtnguyen2/events{/privacy}","received_events_url":"https://api.github.com/users/dtnguyen2/received_events","type":"User","site_admin":false},"created_at":"2016-03-03T11:34:13Z","updated_at":"2016-03-03T11:34:13Z","body":"Sorry for the late response, but it is crunch time for us in another\nproject but I will give you a quick and dirty (i.e. I have not tried to\ncode it but fairly sure that it will work) answer to your question.\n\nThe correct fix is to modify the Fit class, as called in line 161 of your\ncode to pass in the extra statistic etc... but I suspect that is going to\nbe a lot of work.\n\nThe quick and dirty way is the following and I suspect when I get time to\nlook into this I will fix wstat once and for all.  When I added wstat\nrecently, please see the file sherpa/stats/__init__.py.  At the end of the\naforementioned file, one should be able to see the following code:\n\ndef __init__(self, name='wstat'):\n\n        Likelihood.__init__(self, name)\n\n\n    @staticmethod\n\n    def calc_stat(data, model, staterror=None, syserror=None,\n\n                  weight=None, bkg=None):\n\n        if bkg is None or bkg['bkg'] is None:\n\n            raise StatErr('usecstat')\n\n\n        return _statfcts.calc_wstat_stat(data, model,\n\n                                         bkg['data_size'],\n\n                                         bkg['exposure_time'],\n\n                                         bkg['bkg'],\n\n                                         bkg['backscale_ratio'],\n\n                                         truncation_value)\n\n\nAs you can see, since I did not want to modify the other stat and I don't\nneed the bkg data in the stat code I passed in a bunch of extra info, the\nimportant data is the array data_size.  This array should give you the size\nof fermi_data and hess_data as called on line 155.  So this should enable\nus to calculate the different stats for different data sets.  There should\nstill be some work to do but I think it can be done.\n\nOn Wed, Mar 2, 2016 at 2:03 PM, Christoph Deil <notifications@github.com>\nwrote:\n\n> @dtnguyen2 <https://github.com/dtnguyen2> - Thanks for the advice.\n>\n> @adonath <https://github.com/adonath> and I started working on a\n> prototype that does what we want in gammapy/gammapy#479\n> <https://github.com/gammapy/gammapy/pull/479>\n>\n> The main things we're struggling with is:\n>\n>    1.\n>\n>    How to have several datasets with different fit statistics?\n>    We figured out how to have different data and model objects here\n>    <https://github.com/gammapy/gammapy/pull/479/files#diff-f2a0947dbf9e735fa35ccb13368b4d1eR139>,\n>    but we couldn't figure out how to have different stat objects for each with\n>    Sherpa.\n>    2.\n>\n>    We did manage to hack together code where the data and fit statistic\n>    are not part of Sherpa (but computed by the Fermi ScienceTools or some\n>    other software) here:\n>\n>    https://github.com/gammapy/gammapy/pull/479/files#diff-f2a0947dbf9e735fa35ccb13368b4d1eR116\n>    There must be a better way than this hacky double-shim approach, no?\n>\n> —\n> Reply to this email directly or view it on GitHub\n> <https://github.com/sherpa/sherpa/issues/184#issuecomment-191374730>.\n>\n"}}