{"action":"created","issue":{"url":"https://api.github.com/repos/vegandevs/vegan/issues/197","repository_url":"https://api.github.com/repos/vegandevs/vegan","labels_url":"https://api.github.com/repos/vegandevs/vegan/issues/197/labels{/name}","comments_url":"https://api.github.com/repos/vegandevs/vegan/issues/197/comments","events_url":"https://api.github.com/repos/vegandevs/vegan/issues/197/events","html_url":"https://github.com/vegandevs/vegan/issues/197","id":177180479,"number":197,"title":".Call interface for nullmodel simulation","user":{"login":"jarioksa","id":1432933,"avatar_url":"https://avatars.githubusercontent.com/u/1432933?v=3","gravatar_id":"","url":"https://api.github.com/users/jarioksa","html_url":"https://github.com/jarioksa","followers_url":"https://api.github.com/users/jarioksa/followers","following_url":"https://api.github.com/users/jarioksa/following{/other_user}","gists_url":"https://api.github.com/users/jarioksa/gists{/gist_id}","starred_url":"https://api.github.com/users/jarioksa/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/jarioksa/subscriptions","organizations_url":"https://api.github.com/users/jarioksa/orgs","repos_url":"https://api.github.com/users/jarioksa/repos","events_url":"https://api.github.com/users/jarioksa/events{/privacy}","received_events_url":"https://api.github.com/users/jarioksa/received_events","type":"User","site_admin":false},"labels":[{"url":"https://api.github.com/repos/vegandevs/vegan/labels/enchancement","name":"enchancement","color":"207de5"}],"state":"open","locked":false,"assignee":null,"assignees":[],"milestone":null,"comments":3,"created_at":"2016-09-15T13:43:15Z","updated_at":"2016-09-20T08:31:36Z","closed_at":null,"body":"I have started implementing `.Call()` interface to `nullmodel` simulation. Currently we have several compiled **C** functions and we access these with `.C()` interface as set up in `make.commsim()` and used then by `simulate.nullmodel()`. I wanted to see if these functions can be made more efficient by switching from `.C` to `.Call`. \r\n\r\nIn `.C` interface we pass basic **C** structures, typically **R** data as (pointers to) **C** vectors. The passing of arguments and return of results happens via arguments. In `.Call` interface we pass **R** structures and can return **R** structures. I have now implemented two test cases which together cover all variants in `nestedness.c` (i.e., sequential and non-sequential models). For sequential methods I have implemented `.Call` to trial swap, and for non-sequential the `.Call` to quasiswap. The traditional interface to trial swap is:\r\n```R\r\n    out <- array(0L, c(nr, nc, n))\r\n    out[, , 1] <- .C(\"trialswap\", m = x, nr, nc, thin, PACKAGE = \"vegan\")$m\r\n    for (k in seq_len(n - 1)) out[, , k + 1] <- .C(\"trialswap\", \r\n        m = out[, , k], nr, nc, thin, PACKAGE = \"vegan\")$m\r\n    out\r\n```\r\nWe first create an empty 3D array `out`, then fill the first facet (third dimension, `out[,,1]`) with swapped input matrix (`m`), and then we loop over all `n` simulations and fill the next matrix with the swapped previous one with repeated `.C(\"trialswap\", ...)`. The new interface is:\r\n```R\r\n .Call(\"do_tswap\", as.matrix(x), n, thin, PACKAGE = \"vegan\")\r\n```\r\nThat's all. We pass an **R** matrix, create the 3D array within the compiled code, loop over its third facet in compiled code, and finally return the 3D array. The `do_tswap` function will call the very same `trialswap` **C** code as the old interface, but this call happens within compiled code instead of crossing **R**/**C** border.  Using `.Call` interface normally means moving a part of **R** into **C** while keeping the old **C**.\r\n\r\nI would appreciate tests and comments on the new interface. The new interface is now available in branch `nestedness-Calls` of the devel branch of **vegan**.  The new interfaces can be defined in the `nullmodel` command or in `oecosimu`. For instance, the version above can be invoked as:\r\n```R\r\nnm <- nullmodel(BCI, \"dotswap\")\r\noecosimu(BCI, nestedchecker, \"dotswap\", thin=1000)\r\n```\r\nThe new alternative models are `\"dotswap\"` , `\"do2tswap\"` (an alternative implementation of the same) and `\"doquasiswap\"`, and the calls to old interfaces also work (`\"tswap\"`, `\"quasiswap\"`). The old and new implementation should give equal results with the same random number seed so it should be testing of robustness, flexibility and speed.\r\n\r\nAll other nullmodels can be made as a variant of these or even with the same functions with extra argument of `method`, but I want to gain some experience on the code before implementing these for all null models. In particular, I have just started looking at writing code for `.Call()` interface a couple of weeks ago, and I need to check the code thoroughly. Do you have time to look at these, @psolymos ?\r\n\r\nMy first results are that trial swap indeed is much faster with the `.Call()` interface: I have micorbenchmarked 7-fold speed-up for larger data sets (such as `BCI` and `bryceveg`), but the effect varies with platform and data properties. However, I have not found any clear benefits with quasiswap. At best I have had <1% speed-up, but the results are inconsistent and sometimes (though more rarely) the new interface is even slower than the old interface. It seems that quasiswap spends most of its time within compiled code anyway, and mere interface is such a small extra cost that it can hardly be measured. I looked at the working of the **C** code, and it seems that for instance in `BCI` we need usually over 1.5 Million cycles to qusiswap the input matrix into 0/1 matrix. Even the tiny `sipoo` data need 150,000 passes to quasiswap the input matrix. The quasiswap algorithm could be easily parallelized, and that could give a clear advantage, and would probably be easier in **R** than within **C**. Perhaps it makes sense to keep the old interface for computer-heavy non-sequential methods, but implement sequential methods with `.Call`."},"comment":{"url":"https://api.github.com/repos/vegandevs/vegan/issues/comments/248238320","html_url":"https://github.com/vegandevs/vegan/issues/197#issuecomment-248238320","issue_url":"https://api.github.com/repos/vegandevs/vegan/issues/197","id":248238320,"user":{"login":"jarioksa","id":1432933,"avatar_url":"https://avatars.githubusercontent.com/u/1432933?v=3","gravatar_id":"","url":"https://api.github.com/users/jarioksa","html_url":"https://github.com/jarioksa","followers_url":"https://api.github.com/users/jarioksa/followers","following_url":"https://api.github.com/users/jarioksa/following{/other_user}","gists_url":"https://api.github.com/users/jarioksa/gists{/gist_id}","starred_url":"https://api.github.com/users/jarioksa/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/jarioksa/subscriptions","organizations_url":"https://api.github.com/users/jarioksa/orgs","repos_url":"https://api.github.com/users/jarioksa/repos","events_url":"https://api.github.com/users/jarioksa/events{/privacy}","received_events_url":"https://api.github.com/users/jarioksa/received_events","type":"User","site_admin":false},"created_at":"2016-09-20T08:31:36Z","updated_at":"2016-09-20T08:31:36Z","body":"My test with self-built current R-devel on Linux:\r\n![mb](https://cloud.githubusercontent.com/assets/1432933/18662844/a7864e16-7f24-11e6-81c0-9cbbd05e8faa.png)\r\n\r\nThe model parameters were selected so that all examples run in ca. 100 ms in release (2.4-1) version of `nullmodel`. The `do-`model is the one using `.Call` interface, and the plain name is for the old `.C` interface. The parameters for approximately equal running time are highly dependent on the system used and on the data properties (numbers of rows and columns, fills, abundance values). There were also changes in **C** files so that some speed-up to release version came from this. In particular, the RNG generation is faster now. There is a radical improvement in all cases except in non-sequential null models (`quasiswap_count`, `quasiswap`). The `quasiswap` method is particularly slow and cannot be helped with `.Call`.\r\n\r\nHere the command I used in benchmarking with **microbenchmark** package:\r\n```R\r\nmb <- microbenchmark(\r\n    simulate(nullmodel(BCI, \"swap\"), 400, thin=15),\r\n    simulate(nullmodel(BCI, \"doswap\"), 400, thin=15),\r\n    simulate(nullmodel(BCI, \"tswap\"), 400, thin=400),\r\n    simulate(nullmodel(BCI, \"dotswap\"), 400, thin=400),\r\n    simulate(nullmodel(BCI, \"curveball\"), 300, thin=50),\r\n    simulate(nullmodel(BCI, \"docurveball\"), 300, thin=50),\r\n    simulate(nullmodel(BCI, \"swap_count\"), 400),\r\n    simulate(nullmodel(BCI, \"doswap_count\"), 400),\r\n    simulate(nullmodel(BCI, \"abuswap_r\"), 300),\r\n    simulate(nullmodel(BCI, \"doabuswap_r\"), 300),\r\n    simulate(nullmodel(BCI, \"abuswap_c\"), 300),\r\n    simulate(nullmodel(BCI, \"doabuswap_c\"), 300),\r\n    simulate(nullmodel(BCI, \"quasiswap_count\"), 60),\r\n    simulate(nullmodel(BCI, \"doquasiswap_count\"), 60),\r\n    simulate(nullmodel(BCI, \"quasiswap\"), 2),\r\n    simulate(nullmodel(BCI, \"doquasiswap\"), 2))\r\n```"}}
{"actor":{"display_login":"jarioksa"}}