From mboxrd@z Thu Jan 1 00:00:00 1970 Return-Path: Received: (qmail 7164 invoked by alias); 12 May 2013 18:14:48 -0000 Mailing-List: contact gsl-discuss-help@sourceware.org; run by ezmlm Precedence: bulk List-Id: List-Subscribe: List-Archive: List-Post: List-Help: , Sender: gsl-discuss-owner@sourceware.org Received: (qmail 7153 invoked by uid 89); 12 May 2013 18:14:47 -0000 X-Spam-SWARE-Status: No, score=-4.0 required=5.0 tests=BAYES_00,KHOP_THREADED,RP_MATCHES_RCVD,SPF_PASS autolearn=ham version=3.3.1 Received: from ipmx5.colorado.edu (HELO ipmx5.colorado.edu) (128.138.128.235) by sourceware.org (qpsmtpd/0.84/v0.84-167-ge50287c) with ESMTP; Sun, 12 May 2013 18:14:46 +0000 From: Patrick Alken X-IronPort-Anti-Spam-Filtered: true X-IronPort-Anti-Spam-Result: AnwmAJzbj1FH7SMMPGdsb2JhbAANQwoOwGaCbAMBgRkDAQEBATiCVAEBAQQ4QAEQCxgJFg8JAwIBAgExFAYNAQcBAbIgiT2HXY1sXV8Hg1UDiRqUFY4CXoFMJA Received: from c-71-237-35-12.hsd1.co.comcast.net (HELO [192.168.0.33]) ([71.237.35.12]) by smtp.colorado.edu with ESMTP/TLS/DHE-RSA-CAMELLIA256-SHA; 12 May 2013 12:14:44 -0600 Message-ID: <518FDC14.9070600@colorado.edu> Date: Sun, 12 May 2013 18:14:00 -0000 User-Agent: Mozilla/5.0 (X11; Linux x86_64; rv:17.0) Gecko/20130330 Thunderbird/17.0.5 MIME-Version: 1.0 To: Peter Teuben CC: "gsl-discuss@sourceware.org" Subject: Re: Robust linear least squares References: <518D6E3B.8080503@colorado.edu> <518FD7B4.8070100@astro.umd.edu> In-Reply-To: <518FD7B4.8070100@astro.umd.edu> Content-Type: text/plain; charset=ISO-8859-1; format=flowed Content-Transfer-Encoding: 7bit X-SW-Source: 2013-q2/txt/msg00003.txt.bz2 Hi Peter, The most common robust least squares algorithm is called "M-estimation" which is what I've implemented. At each step of the iteration, you calculate the residuals and use a weighting function which is designed to assign large weights to small residuals and small weights to large residuals, so that the large residuals (outliers) contribute less and less to the model at each iteration. At each iteration, you need an estimate of the residual standard deviation, and I am using the Mean-Absolute-Deviation (MAD) of the p largest residuals (where p is the number of model parameters). There are alternatives to computing sigma but the MAD seems to be the most widely used. If you check out the latest repository, have a look at the manual since I've documented everything including a description of the algorithm used. Let me know if you have more questions. Patrick On 05/12/2013 11:56 AM, Peter Teuben wrote: > Patrick > I agree, this is a useful option! > > can you say a little more here how you define robustness. The one I > know takes the quartiles Q1 and Q3 (where Q2 would > be the median), then define D=Q3-Q1 and only uses points between > Q1-1.5*D and Q3+1.5*D to define things like a robust mean and variance. > Why 1.5 I don't know, I guess you could keep that a variable and tinker > with it. > For OLS you can imagine applying this in an iterative way to the Y > values, since formally the errors in X are neglibable compared to those > in Y. I'm saying iterative, since in theory the 2nd iteration could have > rejected points that should have > been part or the "core points". For non-linear fitting this could be a > lot more tricky. > > peter > > > On 05/10/2013 06:01 PM, Patrick Alken wrote: >> Hi all, >> >> I just committed a significant chunk of code related to robust >> linear regression into GSL and mainly wanted to update the other >> developers and any other interested parties. The main idea here is >> that ordinary least squares is very sensitive to data outliers, and >> the robust algorithm tries to identify and downweight outlier points >> so they don't drastically affect the model. I think this is something >> that has been needed in gsl for a while. >> >> I've been developing the code for a while and have been using it >> successfully in my own work, and also validated it pretty extensively >> against the matlab implementation. I still need to make some automated >> tests for it which I should get to next week. >> >> In the meantime, the code is very usable and working so feel free to >> try it out. >> >> Patrick