2012-10-03 34 views
11

我正試圖使用​​apache.commons.math庫來計算將R腳本轉換爲java的庫。我可以用org.apache.commons.math.analysis.interpolation.LoessInterpolator代替R loess嗎?我無法得到相同的結果。R.loess和org.apache.commons.math之間的區別LoessInterpolator

編輯

這裏是創建一個隨機陣列(X,Y)和計算與LoessInterpolator或通過調用R.在結束時黃土的Java程序,結果被打印。

import java.io.*; 
import java.util.Random; 

import org.apache.commons.math.analysis.interpolation.LoessInterpolator; 


public class TestLoess 
    { 
    private String RScript="/usr/local/bin/Rscript"; 
    private static class ConsummeInputStream 
     extends Thread 
     { 
     private InputStream in; 
     ConsummeInputStream(InputStream in) 
      { 
      this.in=in; 
      } 
     @Override 
     public void run() 
      { 
      try 
       { 
       int c; 
       while((c=this.in.read())!=-1) 
        System.err.print((char)c); 
       } 
      catch(IOException err) 
       { 
       err.printStackTrace(); 
       } 
      } 
     } 
    TestLoess() 
     { 

     } 
    private void run() throws Exception 
     { 
     int num=100; 
     Random rand=new Random(0L); 
     double x[]=new double[num]; 
     double y[]=new double[x.length]; 
     for(int i=0;i< x.length;++i) 
      { 
      x[i]=rand.nextDouble()+(i>0?x[i-1]:0); 
      y[i]=Math.sin(i)*100; 
      } 
     LoessInterpolator loessInterpolator=new LoessInterpolator(
      0.75,//bandwidth, 
      2//robustnessIters 

      ); 
     double y2[]=loessInterpolator.smooth(x, y); 

     Process proc=Runtime.getRuntime().exec(
      new String[]{RScript,"-"} 
      ); 
     ConsummeInputStream errIn=new ConsummeInputStream(proc.getErrorStream()); 
     BufferedReader stdin=new BufferedReader(new InputStreamReader(proc.getInputStream())); 
     PrintStream out=new PrintStream(proc.getOutputStream()); 
     errIn.start(); 
     out.print("T<-as.data.frame(matrix(c("); 
     for(int i=0;i< x.length;++i) 
      { 
      if(i>0) out.print(','); 
      out.print(x[i]+","+y[i]); 
      } 
     out.println("),ncol=2,byrow=TRUE))"); 
     out.println("colnames(T)<-c('x','y')"); 
     out.println("T2<-loess(y ~ x, T)"); 
     out.println("write.table(residuals(T2),'',col.names= F,row.names=F,sep='\\t')"); 
     out.flush(); 
     out.close(); 
     double y3[]=new double[x.length]; 
     for(int i=0;i< y3.length;++i) 
      { 
      y3[i]=Double.parseDouble(stdin.readLine()); 
      } 
     System.out.println("X\tY\tY.java\tY.R"); 
     for(int i=0;i< y3.length;++i) 
      { 
      System.out.println(""+x[i]+"\t"+y[i]+"\t"+y2[i]+"\t"+y3[i]); 
      } 
     } 

    public static void main(String[] args) 
     throws Exception 
     { 
     new TestLoess().run(); 
     } 
    } 

彙編& EXEC:

javac -cp commons-math-2.2.jar TestLoess.java && java -cp commons-math-2.2.jar:. TestLoess 

輸出:

X Y Y.java Y.R 
0.730967787376657 0.0 6.624884763714674 -12.5936186703287 
0.9715042030481429 84.14709848078965 6.5263049649584 71.9725380029913 
1.6089216283982513 90.92974268256818 6.269100654071115 79.839773167581 
2.159358633515885 14.112000805986721 6.051308261720918 3.9270340708818 
2.756903911313087 -75.68024953079282 5.818424835586378 -84.9176311089431 
3.090122310789737 -95.89242746631385 5.689740879461759 -104.617807889069 
3.4753114955304554 -27.941549819892586 5.541837854229562 -36.0902352062634 
4.460153035730264 65.6986598718789 5.168028655980764 58.9472823439219 
5.339335553602744 98.93582466233818 4.840314399516663 93.3329030534449 
6.280584733084859 41.21184852417566 4.49531113985498 36.7282165788057 
6.555538699120343 -54.40211108893698 4.395343460231256 -58.5812856445538 
6.68443584999412 -99.99902065507035 4.348559404444451 -104.039069260889 
6.831037507640638 -53.657291800043495 4.295400167908642 -57.5419313320511 
6.854275630124528 42.016703682664094 4.286978656933373 38.1564179414478 
7.401015387322993 99.06073556948704 4.089252482141094 95.7504087842369 
8.365502247999844 65.02878401571168 3.7422883733498726 62.5865641279576 
8.469992934250815 -28.790331666506532 3.704793544880599 -31.145867173504 
9.095139297716374 -96.13974918795569 3.4805388562453574 -98.0047896609079 
9.505935493207435 -75.09872467716761 3.3330472034239405 -76.6664588290508 

y的輸出值明顯不是R和爪哇之間是相同的; Y.R列看起來不錯(它接近原始的Y列)。我應該如何改變這個以獲得Y.java〜Y.R?

+1

可能會在'loess'調用中添加'span = 2/3'?我不知道'loess'中的span是否與'LoessInterpolator'中的'bandwidth'參數相同,但'loess'中'span'的默認值是0.75,並且您將「bandwidth」設置爲2/3 。 –

+0

謝謝大家的回答。我遠離我的工作。我明天會探討你的建議。 – Pierre

+2

快速評論。我只是檢查了R和Java中的lo(w)ess實現是完全相同的。 R中的lowess()和Apache Commons Math中的LoessInterpolator()都指向相同的Cleveland(1979)論文並且具有相同的參數。我從兩個實現中獲得了相同的結果。 – user2065369

回答

5

你需要改變的三個輸入參數的默認值,以使在Java和R版本相同:

  1. Java的LoessInterpolator只做線性局部多項式迴歸,但R支持線性(度= 1) ,二次方程(度數= 2)和一個奇怪的度數= 0選項。因此,您需要在R中指定degree=1以與Java相同。

  2. LoessInterpolator默認迭代次數DEFAULT_ROBUSTNESS_ITERS=2,但R默認iterations=4。所以你需要在R中設置control = loess.control(iterations=X)(X是迭代次數)。

  3. LoessInterpolator默認值DEFAULT_BANDWIDTH=0.3但是R默認span=0.75

3

我不能爲Java實現說話,但lowess有一些它控制擬合的帶寬參數。除非您使用相同的控制參數,否則您應該預期結果會有所不同。每當人們平滑數據時,我的建議是繪製原始數據以及擬合數據,並自行決定哪些控制參數可以在數據保真度和平滑(又稱噪聲消除)之間產生所需的折衷。

1

這裏有兩個問題。首先,如果繪製生成的數據,它看起來幾乎是隨機的,並且由R中的黃土產生的擬合非常差,例如,

plot(T$x, T$y) 
lines(T$s, T2$fitted, col="blue", lwd=3) 

plot of the data generated by the Java code above with a loess fit generated by R

然後在你的[R腳本,你正在寫的殘差不是預言所以在這一行

out.println("write.table(residuals(T2),'', 
    col.names= F,row.names=F,sep='\\t')"); 

您需要更改residuals(T2)predict(T2)例如

out.println("write.table(predict(T2),'', 
    col.names= F,row.names=F,sep='\\t')"); 

因此,在您的代碼示例中,R生成的第一對殘差行看起來非常合適,這是純粹的機會。

對我而言,如果我嘗試使用一些更合適的數據進行擬合,那麼Java和R會返回相似但不完全相同的結果。另外我發現結果更接近,如果我沒有調整默認robustnessIter設置。