On Wed, 2009-07-29 at 20:38 -0400, Skipper Seabold wrote:
> On Wed, Jul 29, 2009 at 8:17 PM, Skipper Seabold<jsseab...@gmail.com> wrote:
> > On Thu, Jul 23, 2009 at 9:50 AM, Skipper Seabold<jsseab...@gmail.com> wrote:
> >> On Thu, Jul 23, 2009 at 7:12 AM, Peter<rpy-l...@maubp.freeserve.co.uk>
> >> wrote:
> >>> On Thu, Jul 23, 2009 at 12:02 PM, Laurent Gautier<lgaut...@gmail.com>
> >>> wrote:
> >>>>>
> >>>>> Since x is a NumPy array these can come in different flavours, and
> >>>>> perhaps something here is breaking rpy2 or R, e.g.
> >>>>
> >>>> The bug report if for rpy-1.0.x, I think.
> >>>
> >>
> >> In [2]: rpy.__version__
> >> Out[2]: '1.4.0.dev7055'
> >>
> >> Is that expected behavior? That is the numpy release that I currently
> >> have installed. Note that I tested this on two machines however. But
> >> yes it's rpy-1.0.3 I believe.
> >>
> >>> Sorry, you are right. I would still try making a copy of the NumPy array
> >>> which should give you a continuous C ordered array which owns its
> >>> own memory.
> >>>
> >>> Peter
> >>>
> >>
> >>
> >> Thanks, but I think I checked this. When I said exactly the same, as
> >> far as I can tell they are exactly the same.
> >>
> >> Say exog is my x that I have imported and x is built as above from rpy
> >> and numpy.
> >>
> >> In [16]: x.flags
> >> Out[16]:
> >> C_CONTIGUOUS : False
> >> F_CONTIGUOUS : True
> >> OWNDATA : False
> >> WRITEABLE : True
> >> ALIGNED : True
> >> UPDATEIFCOPY : False
> >>
> >> In [17]: exog.flags
> >> Out[17]:
> >> C_CONTIGUOUS : False
> >> F_CONTIGUOUS : True
> >> OWNDATA : False
> >> WRITEABLE : True
> >> ALIGNED : True
> >> UPDATEIFCOPY : False
> >>
> >> In [18]: x.flags == exog.flags
> >> Out[18]: True
> >>
> >> In [19]: x.dtype
> >> Out[19]: dtype('float64')
> >>
> >> In [20]: exog.dtype
> >> Out[20]: dtype('float64')
> >>
> >> In [21]: x.shape
> >> Out[21]: (21, 4)
> >>
> >> In [22]: exog.shape
> >> Out[22]: (21, 4)
> >>
> >> I also tried making a copy of exog, so I had a C contiguous array that
> >> owns its own data, but it was the same result. The weird thing is
> >> that using exog works for psi.bisquare and psi.huber, but when I try
> >> psi.hampel I get an error about the conformability of x and the wts,
> >> so it appears to be some sort of a dimensions issue in the iterated
> >> weighted least squares steps for RLM. Anyway, I've update my code to
> >> use the array as built up from rpy for now and it works. When I have
> >> time I will try to dig deeper. I might also try to reinstall rpy
> >> tonight, though again this happened on two machines.
> >>
> >> Cheers,
> >>
> >> Skipper
> >>
> >
> > For posterity's sake, I will post my solution. The problem was in how
> > Hampel's norm works and there are two workarounds. The first is to
> > add an axis to the endogenous variable (as the initial exception said.
> > I just looked in the wrong place.). So in my example above for my
> > variables not from R/Rpy, I would have to add an axis to the response
> > variable so y = y[:,numpy.newaxis] or y = y[:,None]. This needs to be
> > done because Hampel's norm is a hard redescender for large outliers
> > (big ones get a weight of zero rather than just downweighted), so the
> > starting values are chosen differently and it was here that my
> > approach was failing.
> >
> > The second and probably "correct" way to do it (the above lead to some
> > odd behavior in how the results were returned) is to specify how you
> > want the initial values to be calculated. So my command is then
> >
> > results3 = r.rlm(formula, data=frame, psi='psi.hampel', init='lts') #
> > lts = Least Trimmed Squares
> >
> > and everything works as usual.
> >
> > Skipper
> >
>
> I spoke too soon. You still need to add the axis, and the results are
> returned as dicts rather than arrays in my script though nothing else
> is different. The weirdest thing has been the consistency with which
> I can reproduce this...
Did you have better luck with rpy2 ?
------------------------------------------------------------------------------
Let Crystal Reports handle the reporting - Free Crystal Reports 2008 30-Day
trial. Simplify your report design, integration and deployment - and focus on
what you do best, core application coding. Discover what's new with
Crystal Reports now. http://p.sf.net/sfu/bobj-july
_______________________________________________
rpy-list mailing list
rpy-list@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/rpy-list
