venerdì 17 settembre 2010

Come visualizzare il codice delle funzioni di R

Spesso ho avuto la necessità di esaminare il codice di una funzione di R, per studiare come essa funziona.
In genere è sufficiente digitare il nome della funzione, senza parentesi o altri parametri per visualizzarlo sulla console di R, ma non sempre questo è sufficiente.


Consideriamo ad esempio la funzione shapiro.test. Per ottenere il codice con cui viene eseguito questo test è sufficiente digitare:


> shapiro.test
function (x)
{
DNAME <- deparse(substitute(x))
x <- sort(x[complete.cases(x)])
stopifnot(is.numeric(x))
n <- length(x)
if (n < 3 || n > 5000)
stop("sample size must be between 3 and 5000")
rng <- x[n] - x[1L]
if (rng == 0)

... ...
... ...


Ma prendiamo adesso in considerazione la funzione anova. Digitando tale parola riservata di R otteniamo:


> anova
function (object, ...)
UseMethod("anova")


Il codice è apparentemente nascosto. Utilizziamo allora la funzione methods:


> methods(anova)

[1] anova.coxph* anova.coxphlist* anova.Design anova.glm anova.glmlist
[6] anova.glmmPQL* anova.lm anova.loess* anova.loglm* anova.mlm
[11] anova.multinom* anova.negbin* anova.nls* anova.polr* anova.survreg*
[16] anova.survreglist*

Non-visible functions are asterisked


Abbiamo ottenuto un elenco di funzioni, complessivamente raggruppate nella funzione anova. Possiamo visualizzare il codice delle funzioni non asteriscate utilizzando la funzione print():


> print(anova.lm)
function (object, ...)
{
if (length(list(object, ...)) > 1L)
return(anova.lmlist(object, ...))
w <- object$weights
ssr <- sum(if (is.null(w)) object$residuals^2 else w * object$residuals^2)
mss <- sum(if (is.null(w)) object$fitted.values^2 else w *
object$fitted.values^2)
if (ssr < 1e-10 * mss)
warning("ANOVA F-tests on an essentially perfect fit are unreliable")
dfr <- df.residual(object)

... ...
... ...


Se invece vogliamo visualizzare il codice di una funzione asteriscata, con lo stesso metodo otteniamo:


> print(anova.coxph)
Errore in print(anova.coxph) : oggetto "anova.coxph" non trovato


Utilizziamo allora la funzione getAnywhere():


> print(getAnywhere(anova.coxph))
A single object matching ‘anova.coxph’ was found
It was found in the following places
registered S3 method for anova from namespace survival
namespace:survival
with value

function (object, ..., test = "Chisq")
{
dotargs <- list(...)
named <- if (is.null(names(dotargs)))
rep(FALSE, length(dotargs))
else (names(dotargs) != "")
if (any(named))
warning(paste("The following arguments to anova.coxph(..)",
"are invalid and dropped:", paste(deparse(dotargs[named]),
collapse = ", ")))
dotargs <- dotargs[!named]

... ...
... ...