Wolfram的《一种新科学》

Stephen Wolfram

15岁开始发表科技论文，20岁获得博士学位，22岁获美国MacArthur奖（最年轻的获奖者）

1986年，开创Wolfram研究公司，专门开发Mathematica软件

A New Kind of Science

1980年，他的研究兴趣从物理学转向了复杂性科学

1991年，开始写作巨著《A New Kind of Science》

2002年5月，正式发表厚达1000多页的《A New Kind of Science》

该书的分析与证明没用数学公式，而是用简单程序画出的图形

NKS就是研究各类简单程序的理论与方法

复杂性科学的应用

数字艺术、互动游戏

分形艺术

计算机动画

计算机游戏

进化计算、涌现计算（Emergent Computation）、仿生计算、人工生命

遗传算法

PSO粒子群优化

蚁群优化算法

真实世界的机器人工生命

复杂系统的计算机建模

复杂网络

社会系统模拟

涂晓媛的人工鱼

毕业于清华大学自动化系（84级）

涂晓媛获得ACM 图灵奖 - 优秀博士论文

“晓媛的鱼” 轰动欧美计算机届、娱乐届

人工生命与粒子群优化算法

函数优化问题

PSOfunc[x_]:=10*2+Sum[((x[[i]])^2)-10Cos[2 Pi x[[i]]],{i,2}]

Plot3D[PSOfunc[{x,y}],{x,-5.12,5.12},{y,-5.12,5.12}]

ContourPlot[PSOfunc[{x,y}],{x,-5.12,5.12},{y,-5.12,5.12}]

PSOfunc[{0,0}]

Boid模型与粒子群优化算法 （Particle Swarm Optimization, PSO）

1995 年，James Kennedy和Russell C.Eberhart受到Boid模型启发，提出PSO算法

基本思想：让一群Boid在问题空间中觅食，函数值越小，表示食物越多

当一只Boid遇到一个食物富足的地方的时候，它会以一定的概率停在那里，并招呼其它Boid过来

Boid又会以一定的概率四处乱飞去觅食

用PSO求解函数优化问题

建模思路：用一个列表存储所有Boid的位置和速度，用gBest存储最多的食物点，用lBest存储每个Boid遇到的最多食物点

初始化

BoidsNum=50;
Inf=10000;
PSOStates=Table[{{RandomReal[{-5.12,5.12}],RandomReal[{-5.12,5.12}]},{0,0}},{i,BoidsNum}];
gBest=Inf;
lBests=Table[Inf,{i,BoidsNum}];
XgBest={-1,1};
XlBests=Table[{10,10},{BoidsNum}];
w=0.05;
C1=0.05;
C2=0.1;

每一步

PSOOnestep[lastState_]:=Module[{i,minFitness,NewV=Table[0,{BoidsNum}],NewX=Table[0,{BoidsNum}],out},
For[i=1,i≤BoidsNum,i++,
If[lBests[[i]]>PSOfunc[lastState[[i,1]]],
lBests[[i]]=PSOfunc[lastState[[i,1]]];
XlBests[[i]]=lastState[[i,1]];
];
];
minFitness=Min[lBests];
If[minFitness<gBest,
gBest=minFitness;
XgBest=lastState[[Ordering[lBests,1][[1]],1]];
];

For[i=1,i≤BoidsNum,i++,
NewV[[i]]=w lastState[[i,2]]+C1 RandomReal[{0,1}] (XlBests[[i]]-lastState[[i,1]])+C2 RandomReal[{0,1}] (XgBest-lastState[[i,1]]);
(*三种力：惯性、趋向我找到的最优掉、趋向整体最优点*)
NewX[[i]]=Mod[lastState[[i,1]]+NewV[[i]]+3*5.12,2*5.12]-5.12;
];
Print[{XgBest,gBest}];
Return[Transpose[{NewX,NewV}]];
];
DrawBoid[state_]:=Arrow[{state[[1]],state[[1]]+state[[2]]}]

做成动画

Module[{last=1},Animate[If[last≠t,last=t;PSOStates=PSOOnestep[PSOStates]];Show[ContourPlot[PSOfunc[{x,y}],{x,-5.12,5.12},{y,-5.12,5.12},Contours->{Automatic,5}],Graphics[DrawBoid[#]&/@PSOStates,PlotRange ->  {{-5.12, 5.12}, {-5.12, 5.12}}]],{t,1,Infinity},AnimationRunning->False,AnimationRate->0.000005]]

复杂网络

网络无处不在

各类网络有一些相似的性质：小世界、无尺度

网络小世界

一个常见的现象：

在一次聚会中，你和一位初次见面的朋友闲聊起来，没多久，你们惊喜地发现，你们都认识某某人，于是感叹： “这个世界真小啊！”

Milgram的“小世界”试验

1961，耶鲁大学的社会心理学家Stanley Milgram做了一个有趣的试验

试验发现，每两个人平均间隔5个人就能发生联系，这个现象被叫做“六度分离”

小世界现象普遍存在于复杂网络中

网络任意两个节点之间的平均距离很小

网络具有聚集效应：你的朋友们也彼此相互熟识

集智俱乐部网络社区www.swarmagents.com

集智俱乐部社区网络分析

从数据库提取网络

file1=#[[1]]&/@Import["dbdata1.dat"];
file2=StringSplit[Import["dbdata2.dat","Text"]];
file3=#[[1]]&/@Import["dbdata3.dat"];
OriDT1={file1,file2,file3};

Poss=Position[OriDT1[[3]],#]&/@Union[DeleteCases[OriDT1[[3]],0.]];
Authors=StringReplace[#,"("~~__->""]&/@OriDT1[[2]];
AuthorGroup=Table[Authors[[#[[1]]]]&/@Poss[[i]],{i,1,Length[Poss]}];
net=Flatten[DeleteCases[ReplaceAll[#,List->Rule]&/@Partition[Union[#],2,1]&/@AuthorGroup,{}]];
netList=net/.Rule->List;
Nodes=Union[Transpose[netList][[1]],Transpose[netList][[2]]];
Length[Nodes]
Length[netList]

网络动态

Manipulate[GraphPlot3D[Take[net,i],MultiedgeStyle->False,VertexLabeling->b,Method->MS],{i,1,Length[net],1},{b,{True,False}},{MS,{"SpringEmbedding","SpringElectricalEmbedding","HighDimensionalEmbedding","RandomEmbedding"}}]

计算平均距离

AverageDistance[netList]

计算聚类系数

N[ClusteringCoefficiency[netList]]

N[ClusteringCoefficiency[Table[{RandomInteger[{0,567}],RandomInteger[{0,567}]},{2213}]]]

集智俱乐部社区网络是一个 “小世界” 网络

图形化小世界现象

DynamicModule[{rep,lastName="",lastType=0},Manipulate[If[Name=="",GraphPlot[net,MultiedgeStyle->False,Method->MS,PlotStyle->Directive[Black],VertexRenderingFunction->Function[{p,l},{Red,Point[p]}]],If[Name≠lastName||lastType≠ShowType,lastType=ShowType;lastName=Name;rep=MyNetworkPlot[netList,Nodes,Name,ShowType]];Graphics[Join[{Black},Line/@DrawNetworks[netList,rep],{Black},{Text[Name,{0,0.4},Background->White]},{Red},Disk[#,0.1]&/@rep[[2]]]]],{Name,Join[{""},Reverse[Nodes]]},{{ShowType,0},{0,1}},{MS,{"SpringElectricalEmbedding","SpringEmbedding","HighDimensionalEmbedding","RandomEmbedding"}}]]

计算度分布

nLinks=Table[{Nodes[[i]],Select[netList,#[[1]]==Nodes[[i]]||#[[2]]==Nodes[[i]]&]},{i,1,Length[Nodes]}];
degrees=Length[#[[2]]]&/@nLinks;
Tally[degrees];
ListPlot[Tally[degrees],AxesLabel->{"连接数","节点数"}]
fit=Fit[DeleteCases[{Log[#[[1]]],Log[#[[2]]]}&/@Tally[degrees],{___,0}],{1,x},x];
Show[ListLogLogPlot[Tally[degrees],AxesLabel->{"Log[连接数]","Log[节点数]"}],Plot[fit,{x,0,1000},PlotStyle->{Red}]]

{#[[1]],Length[#[[2]]]}&/@Sort[nLinks,Length[#1[[2]]]>Length[#2[[2]]]&]//ColumnForm

结束语

简单的程序孕育着复杂的行为

通过计算机模拟，我们可以更好地理解复杂系统

“给我一个支点，我能撬起整个地球” ——阿基米德

“给我一台计算机，我能创造一个世界” ——张江

Mathematica是一个强大、便利的探索平台

初始化

定义系统变量

AppendTo[$Path,NotebookDirectory[]<>"Data"];
AppendTo[$Path,NotebookDirectory[]];

网络的平均路径函数

AverageDistance[network_]:=Module[{SNodes,TheMatrix,AllMatrix,DisMatrix,i,j,k,temp1,temp2},
SNodes=Union[First[Transpose[network]],First[Rest[Transpose[network]]]];
DisMatrix=Table[0,{i,Length[SNodes]},{j,Length[SNodes]}];
TheMatrix=Table[If[Cases[network,{SNodes[[i]],SNodes[[j]]}]!={}||Cases[network,{SNodes[[j]],SNodes[[i]]}]≠{},1,0],{i,1,Length[SNodes]},{j,1,Length[SNodes]}];
AllMatrix=TheMatrix;
For[i=1,i<Length[SNodes],i++,
temp1=0;
temp2=1;
For[j=1,j<=Length[SNodes],j++,
temp2=0;
For[k=1,k<=Length[SNodes],k++,
If[DisMatrix[[j,k]]==0&&AllMatrix[[j,k]]≠0&&j≠k,DisMatrix[[j,k]]=i;temp2=1];
];
temp1+=temp2;
];
If[temp1==0,Break[]];
AllMatrix=AllMatrix .TheMatrix;
];
(*Print[DisMatrix//MatrixForm];*)
Return[N[Mean[Mean/@DisMatrix]]];
];

聚类系数

ClusteringCoefficiency[network_]:=Module[{Neighbors,SNodes,TheMatrix,Clusters,Count=0},
SNodes=Union[First[Transpose[network]],First[Rest[Transpose[network]]]];
Neighbors=Table[DeleteCases[Union[Flatten[Select[network,#[[1]]==SNodes[[i]]||#[[2]]==SNodes[[i]]&]]],SNodes[[i]]],{i,1,Length[SNodes]}];
(*Print[Neighbors];*)
Clusters=(If[Length[#]≤1,0,(Count=0;For[i=1,i≤Length[#]-1,i++,For[j=i+1,j≤Length[#],j++,If[Union[Cases[network,{#[[i]],#[[j]]}],Cases[network,{#[[j]],#[[i]]}]]≠{},Count++;];];];Count
)/(Length[#](Length[#]-1))])&/@Neighbors;
Print[Mean[Clusters]];
Return[N[Mean[Clusters]]];
]

自定义网络图

MyNetworkPlot[Network_,Nodes_,CenterNode_,tp_]:=Module[{lsts=Table[{0,0},{Length[Nodes]}],Visited={},temp,temp1,temp2,temp3,rnd},temp=NestList[(Flatten[Map[DeleteCases[Union[Flatten[Transpose[Cases[Network,{#,_}|{_,#}]]]],#]&,#]])&,{CenterNode},3];
temp2=DeleteCases[Table[temp1=Complement[temp[[i]],Visited];Visited=Union[Visited,temp[[i]]];temp1,{i,1,Length[temp]}],{}];
(*temp3=Transpose[Nodes][[1,Position[Transpose[Nodes][[1]],#][[1,1]]]]&/@Flatten[temp2];*)
If[tp==0,
lsts=Flatten[Table[rnd=RandomReal[{0,1}];Table[{2(i-1) Cos[2 (j+rnd) Pi/Length[temp2[[i]]]],2 (i-1) Sin[2 (j+rnd) Pi/Length[temp2[[i]]]]},{j,Length[temp2[[i]]]}],{i,Length[temp2]}],1];
,
lsts=Flatten[Table[Table[{500(i-1),2j},{j,Length[temp2[[i]]]}],{i,Length[temp2]}],1];
];
Return[{Flatten[temp2],lsts}];
]
DrawNetworks[Network_,NodesRep_]:=Module[{partNet},partNet=DeleteCases[Table[If[Cases[NodesRep[[1]],Network[[i,1]]]≠{}&&Cases[NodesRep[[1]],Network[[i,2]]]≠{},Network[[i]]],{i,Length[Network]}],Null];(partNet/.({x_,y_}:>{Position[NodesRep[[1]],x][[1,1]],Position[NodesRep[[1]],y][[1,1]]}))/.({x_,y_}:>{NodesRep[[2,x]],NodesRep[[2,y]]})]