Thu 07 Nov 2013 05:41:05 PM UTC, comment #30:
I checked in a final changeset to have the moment function default to central moments instead of raw moments (http://hg.savannah.gnu.org/hgweb/octave/rev/58b39152b0f6). All three changes to moment, skewness, and kurtosis will be part of version 3.8.0. Closing report.
|
Sat 19 Oct 2013 02:39:09 PM UTC, comment #29:
I checked in your changeset here (http://hg.savannah.gnu.org/hgweb/octave/rev/f0e777cf348f). I also announced in the NEWS file that the definition of the kurtosis function has changed from "excess kurtosis" to just kurtosis. The NEWS file is a ChangeLog for user-visible changes that are important or that might change behavior from previous versions.
|
Sat 19 Oct 2013 10:13:20 AM UTC, comment #28:
Here is a new version of the patch :
- the tests have been modified as discussed,
- the documentation has been updated.
(file #29418)
|
Sat 19 Oct 2013 06:27:12 AM UTC, comment #27:
For compatibility with MAtlab's kurtosis function, I have changed the definition kurtosis (x).
The new value of kurtosis (x) is equal to the former one ("excess kurtosis") + 3.
Should I mention that in the documentation ? or somewhere else ?
|
Sat 19 Oct 2013 06:14:01 AM UTC, comment #26:
Ok, I'll raise te tolerance.
|
Fri 18 Oct 2013 08:23:55 PM UTC, comment #25:
That is strange, but I'm inclined to accept the results of the arbitrary precision arithmetic. Raising the tolerance is fine.
|
Fri 18 Oct 2013 08:00:59 PM UTC, comment #24:
This is weird, because I get
5.437731792528892
on my computer with both R 3.0.0, Matlab R2012a and my latest Octave build. On the other hand, using mpmath with 30 digits of precision I get the following interval:
mpi('5.43773179252889008331093791989766', '5.43773179252889008331093791992922')
which means that the value that you report is the correct one. I will change the current reference values for more precise ones:
2.9786509002956195 (sample kurtosis)
5.4377317925288901 (with bias corrections)
but then I need to raise the tolerance to 8 * eps on my computer to get the tests to pass. Is that OK ?
Here is the Python code that I use:
|
Fri 18 Oct 2013 05:41:02 PM UTC, comment #23:
I ran 'test kurtosis' and I get a failing test
How did you determine the expected value for this test? Did you use an arbitrary precision calculator? The results for Octave show a difference of '2' in the last significant digit.
|
Fri 18 Oct 2013 05:20:27 PM UTC, comment #22:
Here is a new version of the patch with the two modifications that you suggest. I will fix the documentation during the week-end.
(file #29407)
|
Fri 18 Oct 2013 03:20:14 PM UTC, comment #21:
The patch looks good. I might change
to
Since it "looks" more like a straightforward equation. You don't have to worry about the values not in idx because they have been assigned the value of NaN, and any further math operation with a NaN value produces a NaN.
Also, this was my mistake in one of the tests for division-by-zero. I save and restore the warning state but don't actually force the "divide-by-zero" state to on. I've corrected that in skewness so that the test reads
|
Fri 18 Oct 2013 11:29:43 AM UTC, comment #20:
Here is a first version of my patch for kurtosis.
It's not ready to be pushed yet, the documentation needs to be rewritten, but the code and the tests are ready for review.
@++
(file #29401)
|
Mon 14 Oct 2013 07:02:25 AM UTC, comment #19:
Hello Rik,
I will do my best to provide a similar patch for kurtosis before Friday.
@++
Julien
|
Mon 14 Oct 2013 01:23:53 AM UTC, comment #18:
Julien,
We are getting ready for the 3.8 release. If you have any further patches for kurtosis or moment can I get them within the next week?
Thanks,
Rik
|
Mon 07 Oct 2013 08:14:02 PM UTC, comment #17:
I added the fix to the comment here (http://hg.savannah.gnu.org/hgweb/octave/rev/a53a97faf61a).
|
Mon 07 Oct 2013 07:12:19 PM UTC, comment #16:
Oooops, here is the patch, sorry about that.
(file #29309)
|
Mon 07 Oct 2013 07:11:45 PM UTC, comment #15:
Another tiny patch on skewness.m (Fix comment).
|
Mon 07 Oct 2013 02:52:31 PM UTC, comment #14:
Julien,
I added your changes here (http://hg.savannah.gnu.org/hgweb/octave/rev/8334144458a4).
I thought there was a potential for a "division-by-zero" error which is why I had coded skewness the way it was a long time ago. I added a %!test to verify that if the code is ever changed again we don't re-introduce the same error.
|
Mon 07 Oct 2013 01:28:00 PM UTC, comment #13:
Rik,
Please find attached three patches for skewness.m :
a) Your coding style fixes make some test fail on my computer (that's because you replaced some explicite tolerance values by eps). I simply replace eps by 2 * eps.
b) Because of some changes you made in input argument checking, weird (non scalar) values of flag were able to go though undetected. I fixed that in a manner, I hope, you will find in agreement with Octave core coding conventions.
c) Get rid of a "division by zero" warning that I hadn't noticed before.
@++
Julien
(file #29306, file #29307, file #29308)
|
Mon 07 Oct 2013 02:58:03 AM UTC, comment #12:
Julien,
I reviewed your patch and committed it here (http://hg.savannah.gnu.org/hgweb/octave/rev/ab75e72c5b36). I followed it up with a changeset which modified the coding style to reflect core Octave (http://hg.savannah.gnu.org/hgweb/octave/rev/7004c733412f). I'm attaching the diff file to this bug report which I annotated with the reason for each change. I think that will make it useful for learning. It should be straightforward to do kurtosis now, although I would use the correction factor as shown by the Mathworks in their documentation (http://www.mathworks.com/help/stats/kurtosis.html).
Cheers,
Rik
(file #29303)
|
Sun 06 Oct 2013 10:35:45 AM UTC, comment #11:
Hello Rik,
I can make a similar patch for kurtosis, if this one is OK.
Please tell me whenever you have time to look at it.
@++
Julien
|
Wed 02 Oct 2013 09:19:03 PM UTC, comment #10:
Here is another version of the same patch, with a better documentation (the TeX output looks ok now).
(file #29288)
|
Wed 02 Oct 2013 06:26:26 PM UTC, comment #9:
Ok, I'm trying now to improve the TeX documentation.
|
Wed 02 Oct 2013 05:20:14 PM UTC, comment #8:
Julien,
The TeX part shows up in:
octave/doc/interpreter/octave.pdf
This is the Octave Manual in PDF form.
|
Wed 02 Oct 2013 11:11:04 AM UTC, comment #7:
Ok, here is a first version of my patch for skewness.m, for review.
It solves several problems:
1) Add the missing 'flag' argument for Matlab compatibility
2) Return NaN when the standard deviation is zero (also for Matlab compatibility)
3) Deal when the type of the input argument (single/double) in a more transparent manner, in order to make it easir to solve bug #32168 (related to dataframes).
I have updated the documentation, but I probably need some help to make it right, since this is my first attempt at writing texinfo doc. In particular, I have tried to update the tex-specific part of the doc, but I don't known where I should check that the result is correct...
@++
Julien
(file #29282)
|
Tue 01 Oct 2013 10:53:15 PM UTC, comment #6:
Execellent, I'll review it whenever it's ready.
|
Tue 01 Oct 2013 09:25:21 PM UTC, comment #5:
Hello Rik,
I have a patch almost ready for skewness.m that fixes the problem discussed here, fixes another compatibility issue (skewness should return NaN when the standard deviation is zero) and should make it easy to close bug #32168 at the same time.
I just need to update the documentation, I'll upload the patch in a few days at most.
@++
Julien
|
Thu 26 Sep 2013 10:34:06 PM UTC, comment #4:
I had simply forgotten about this. I'll see if I can fix this in within the next week or so.
|
Thu 26 Sep 2013 07:54:33 PM UTC, comment #3:
I checked now with version 3.7.7+ and the problem is still present. The status is "Need Info". Is there more information that I can provide?
|
Sun 01 Jul 2012 07:30:00 PM UTC, comment #2:
Here are the results with bias-correction:
>> x = [1:5].^2;
>> skewness (x, 0)
ans =
0.6996
>> kurtosis (x, 0)
ans =
2.3577
Thanks and best regards,
Alexander
|
Fri 29 Jun 2012 06:13:30 PM UTC, comment #1:
I swear I updated this a few years ago. Maybe Matlab has been slowly changing their defaults.
On reading the documentation I see that skewness and kurtosis have a flag which decides whether bias correction is applied. Rather than just changing Octave to use the current default, it would be better to add support for the bias flag and then set the flag to the default value used by Matlab.
Could you run the following tests to see how Matlab is applying bias correction?
|
Sun 24 Jun 2012 08:46:06 PM UTC, original submission:
In matlab (R2010a), I get:
>> x = [1:5].^2;
>> moment(x,3)
ans =
303.6000
>> skewness(x)
ans =
0.4693
% this is (std(x,1)^(-3)) * mean ((x - mean(x)).^3)
>> kurtosis(x)
ans =
1.8394
while in octave 3.6.0:
moment(x,3)
ans = 4103
In matlab, moments are central while octave uses raw moment per default (central moments can be activated with moment(x,3,'c'))
skewness(x)
ans = 0.33580
This is (std(x,0)^(-3)) * mean ((x - mean(x)).^3). The normalization of std is thus different.
kurtosis(x)
ans = -1.8228
Again a different normalization of std, and a minus 3 in octave (excess kurtosis).
Concerning the normalization of std (1/n versus 1/(n-1)), I have seen at least a third definition on the web.
(http://www.itl.nist.gov/div898/handbook/eda/section3/eda35b.htm). What a mess!
I would propose that we follow the matlab convention here.
Attached are some patches, with the necessary changes in the calculations to the functions (but help text and test are not updated yet, ... bed time now).
|