[Rd] [External] Re: svd() results should have a class

Lenth, Russell V ru@@e||-|enth @end|ng |rom u|ow@@edu
Fri Jun 24 16:41:08 CEST 2022


Interestingly, for testing purposes, I created an SVD() function exactly like the one you show, even the same name.
What I think I'll do in the estimability package is just quietly add (but not export) a nonest.basis.svd function. Its main use, actually, will be that the nonest.basis.matrix method will be modified to pass most of its work to nonest.basis.svd(svd(x)) since, as it turns out, svd() is the simplest and most efficient way to obtain the desired result -- better than what I have now. I may also add a .list method that checks to see if the right named components of the right classes for svd are available, then passes to the .svd method.

So what I am saying is that there is no urgency to making the requested change.

If you do decide to add a class to svd() values (and La.svd()??), I defer to you experts on whether to inherit from "list". 



-----Original Message-----
From: John Fox <jfox using mcmaster.ca> 
Sent: Friday, June 24, 2022 8:56 AM
To: Lenth, Russell V <russell-lenth using uiowa.edu>
Cc: r-devel using r-project.org; Duncan Murdoch <murdoch.duncan using gmail.com>; Robert Harlow <rharlow86 using gmail.com>
Subject: Re: [Rd] [External] Re: svd() results should have a class

Dear Russ,

A couple of thoughts:

First, it's not clear to me where the svd objects are coming from -- that is, are they generated in the estimability package or coming, e.g., from model objects? If the former, then you could do something like this, creating your own svd class:

 > SVD <- function(x, ...){
+   s <- svd(x, ...)
+   class(s) <- "svd" # or c("svd", class(s))
+   s
+ }

And then, e.g.,

 > kappa.svd <- function(z, tol=sqrt(.Machine$double.eps), ...){
+   d <- z$d
+   d <- d[abs(d) > tol]
+   d[1]/d[length(d)]
+ }

 > m <- lm(Employed ~ ., data=longley)
 > kappa(m, exact=TRUE)
[1] 23845862
 > kappa(SVD(model.matrix(m)))
[1] 23845862

(I suspect you would have thought of this if it solved your problem.)

Second, I think that it's a good idea more generally to have base::svd() return objects of class "svd", and I doubt whether it really matters to append "list" to the class vector. It's not impossible, of course, that changing the class to "svd" from "list" would break something, but it's also hard to imagine a .list method that would expect an svd object.


On 2022-06-23 8:45 p.m., Lenth, Russell V wrote:
> Duncan et al.,
> The real example would be the estimability package, in which there is a generic function nonest.basis and methods for classes lm, matrix, and qr. I wanted to add a method for class svd except there is no such class.
> I agree BTW that the safe thing to do would be to have the returned object be of class c("svd", "list").
> Best,
> Russ
> -----Original Message-----
> From: Duncan Murdoch <murdoch.duncan using gmail.com>
> Sent: Thursday, June 23, 2022 7:23 PM
> To: Lenth, Russell V <russell-lenth using uiowa.edu>; Robert Harlow 
> <rharlow86 using gmail.com>
> Cc: r-devel using r-project.org
> Subject: Re: [Rd] [External] Re: svd() results should have a class
> On 23/06/2022 8:07 p.m., Lenth, Russell V wrote:
>> Bob,
>> I'm not talking about using svd as a generic method. I am talking 
>> about a method FOR svd results, e.g. an S3 method like foo.svd(), for 
>> which there already exist other methods, say foo.default and foo.qr.
>> Currently if I wanted to do
>>       svdobj <- svd(x)
>>       foo(svdobj)
>> it would not dispatch correctly because there is no svd class. Instead, it would be handled by foo.list if it exists, and it is certainly not clear that foo.list would do the right thing.
> I think this would be more convincing if you gave a real example.  You 
> can see the existing methods for the "qr" and "eigen" classes using
>     methods(class = "qr")
>     methods(class = "eigen")
> When I do it, I think I'm only seeing methods from base packages, and they are:
>     kappa for qr
>     solve for qr
>     print for eigen
> I think the default print method works fine for svd() results.
> kappa.svd would probably make sense, but would need some thought:
> calculating kappa on a matrix and using the help page for kappa to naively calculate it from the svd() of that matrix give different results:
>     > set.seed(123)
>     > X <- matrix(rnorm(25), 5,5)
>     > kappa(X)
>     [1] 90.71283
>     > s <- svd(X)
>     > s$d[1]/s$d[length(s$d)]
>     [1] 62.68048
> Duncan Murdoch
>> Russ
>> Sent from my iPad
>> On Jun 23, 2022, at 6:53 PM, Robert Harlow <rharlow86 using gmail.com> wrote:
>> Don't have a view on whether it makes sense in base R or not, but WRE section 7.1 may be helpful to you: https://cran.r-project.org/doc/manuals/R-exts.html#Adding-new-generics.
>> It's not uncommon for packages to want to make base methods generic and the above link provides advice on how to do so.
>> Bob
>> On Thu, Jun 23, 2022 at 12:07 PM Lenth, Russell V <russell-lenth using uiowa.edu<mailto:russell-lenth using uiowa.edu>> wrote:
>> Dear R-Devel,
>> I noticed that if we run base::svd(x), we obtain an object of class "list". Shouldn't there be an "svd" class, in case someone (e.g., me) wants to write methods for singular value decompositions? Note that other matrix-decomposition routines like qr() and eigen() each return objects having those names.
>> Thanks
>> Russ Lenth
>> Russell-lenth using uiowa.edu<mailto:Russell-lenth using uiowa.edu>
>>           [[alternative HTML version deleted]]
>> ______________________________________________
>> R-devel using r-project.org<mailto:R-devel using r-project.org> mailing list 
>> https://stat.ethz.ch/mailman/listinfo/r-devel
>> 	[[alternative HTML version deleted]]
>> ______________________________________________
>> R-devel using r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-devel
> ______________________________________________
> R-devel using r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel
John Fox, Professor Emeritus
McMaster University
Hamilton, Ontario, Canada
web: https://socialsciences.mcmaster.ca/jfox/

More information about the R-devel mailing list